This work addresses the challenges of distributing large physics-based simulations often encountered in the visual effects industry. These simulations, based on partial differential equations, model complex phenomena such as free surface liquids, flames, and explosions, and are characterized by domains whose shapes and topologies evolve rapidly. In this context, we propose a novel partitioning algorithm employing optimal transport —which produces a power diagram—and designed to handle a vast variety of simulation domain shapes undergoing rapid changes over time. Our Power partitioner ensures an even distribution of computational tasks, reduces inter-node data exchange, and maintains temporal consistency, all while being intuitive and artist-friendly. To quantify partitioning quality we introduce two metrics, the surface index and the temporal consistency index, which we leverage in a range of comparisons on real-world film production data, showing that our method outperforms the state of the art in a majority of cases.

We propose a parameter optimization method for achieving static equilibrium of discrete elastic rods. Our method simultaneously optimizes material stiffness and rest shape parameters under box constraints to exactly enforce zero net forces while avoiding stability issues and violations of physical laws. For efficiency, we split our constrained optimization problem into primal and dual subproblems via the augmented Lagrangian method, while handling the dual maximization subproblem via simple vector updates. To efficiently solve the box-constrained primal minimization subproblem, we propose a new active-set Cholesky preconditioner for variants of conjugate gradient solvers with active sets. Our method surpasses prior work in generality, robustness, and speed.

Abstract We propose a new rest shape optimization framework to achieve sag‐free simulations of discrete elastic rods. To optimize rest shape parameters, we formulate a minimization problem based on the kinetic energy with a regularizer while imposing box constraints on these parameters to ensure the system's stability. Our method solves the resulting constrained minimization problem via the Gauss‐Newton algorithm augmented with penalty methods. We demonstrate that the optimized rest shape parameters enable discrete elastic rods to achieve static equilibrium for a wide range of strand geometries and material parameters.

Abstract We introduce a novel unified mixture‐based method for simulating underwater bubbles across a range of bubble scales. Our approach represents bubbles as a set of Lagrangian particles that are coupled with the surrounding Eulerian water volume. When bubble particles are sparsely distributed, each particle, typically smaller than the liquid grid voxel size, corresponds to an individual spherical bubble. As the sub‐grid particles increase in local density our model smoothly aggregates them, ultimately forming connected, fully aerated volumetric regions that are properly resolved by the Eulerian grid. We complement our scheme with a continuous surface tension model, defined via the gradient of the bubbles' local volume fractions, which works seamlessly across this scale transition. Our unified representation allows us to capture a wide range of effects across different scales—from tiny dispersed sub‐grid air pockets to fully Eulerian two‐phase interfacial flows.

We propose a new barrier-based box-constrained convex QP solver based on a primal-dual interior point method to efficiently solve large-scale pressure Poisson problems with non-negative pressure constraints, which commonly arise in liquid animation. The performance of prior active-set-based approaches is limited by the need to repeatedly update the active set. Our solver eliminates this issue by entirely avoiding the use of an active set, which in turn makes the inner problems of our Newton iteration process fully unconstrained. For efficiency, exploiting the solution uniqueness of convex QPs and the fact that the pressure constraints are simple box constraints, we aggressively update solution candidates without performing any step selection procedure (such as line search) and instead directly clamp candidates back to the bounds wherever constraint violations occur. Additionally, to accelerate the inner linear solves, we present a topology-aware geometry-inspired aggregation algebraic multigrid preconditioner and describe in detail several key performance optimizations that we incorporate. We demonstrate the efficacy of our solver in various practical scenarios and show that it often surpasses various alternatives in terms of speed and memory usage.

We present projected walk on spheres (PWoS), a novel pointwise and discretization-free Monte Carlo solver for surface PDEs with Dirichlet boundaries, as a generalization of the walk on spheres method (WoS) [Muller 1956; Sawhney and Crane 2020]. We adapt the recursive relationship of WoS designed for PDEs in volumetric domains to a volumetric neighborhood around the surface, and at the end of each recursion step, we project the sample point on the sphere back to the surface. We motivate this simple modification to WoS with the theory of the closest point extension used in the closest point method. To define the valid volumetric neighborhood domain for PWoS, we develop strategies to estimate the local feature size of the surface and to compute the distance to the Dirichlet boundaries on the surface extended in their normal directions. We also design a mean value filtering method for PWoS to improve the method's efficiency when the surface is represented as a polygonal mesh or a point cloud. Finally, we study the convergence of PWoS and demonstrate its application to graphics tasks, including diffusion curves, geodesic distance computation, and wave propagation animation. We show that our method works with various types of surfaces, including a surface of mixed codimension.

Abstract An ongoing challenge in fluid animation is the faithful preservation of vortical details, which impacts the visual depiction of flows. We propose the Impulse Particle‐In‐Cell (IPIC) method, a novel extension of the popular Affine Particle‐In‐Cell (APIC) method that makes use of the impulse gauge formulation of the fluid equations. Our approach performs a coupled advection‐stretching during particle‐based advection to better preserve circulation and vortical details. The associated algorithmic changes are simple and straightforward to implement, and our results demonstrate that the proposed method is able to achieve more energetic and visually appealing smoke and liquid flows than APIC.

Many geometry processing techniques require the solution of partial differential equations (PDEs) on manifolds embedded in ℝ 2 or ℝ 3 , such as curves or surfaces. Such manifold PDEs often involve boundary conditions (e.g., Dirichlet or Neumann) prescribed at points or curves on the manifold’s interior or along the geometric (exterior) boundary of an open manifold. However, input manifolds can take many forms (e.g., triangle meshes, parametrizations, point clouds, implicit functions, etc.). Typically, one must generate a mesh to apply finite element-type techniques or derive specialized discretization procedures for each distinct manifold representation. We propose instead to address such problems in a unified manner through a novel extension of the closest point method (CPM) to handle interior boundary conditions. CPM solves the manifold PDE by solving a volumetric PDE defined over the Cartesian embedding space containing the manifold and requires only a closest point representation of the manifold. Hence, CPM supports objects that are open or closed, orientable or not, and of any codimension. To enable support for interior boundary conditions, we derive a method that implicitly partitions the embedding space across interior boundaries. CPM’s finite difference and interpolation stencils are adapted to respect this partition while preserving second-order accuracy. Additionally, we develop an efficient sparse-grid implementation and numerical solver that can scale to tens of millions of degrees of freedom, allowing PDEs to be solved on more complex manifolds. We demonstrate our method’s convergence behavior on selected model PDEs and explore several geometry processing problems: diffusion curves on surfaces, geodesic distance, tangent vector field design, harmonic map construction, and reaction-diffusion textures. Our proposed approach thus offers a powerful and flexible new tool for a range of geometry processing tasks on general manifold representations.

Signed distance fields (SDFs) are a widely used implicit surface representation, with broad applications in computer graphics, computer vision, and applied mathematics. To reconstruct an explicit triangle mesh surface corresponding to an SDF, traditional isosurfacing methods, such as Marching Cubes and and its variants, are typically used. However, these methods overlook fundamental properties of SDFs, resulting in reconstructions that exhibit severe oversmoothing and feature loss. To address this shortcoming, we propose a novel method based on a key insight: each SDF sample corresponds to a spherical region that must lie fully inside or outside the surface, depending on its sign, and that must be tangent to the surface at some point. Leveraging this understanding, we formulate an energy that gauges the degree of violation of tangency constraints by a proposed surface. We then employ a gradient flow that minimizes our energy, starting from an initial triangle mesh that encapsulates the surface. This algorithm yields superior reconstructions to previous methods, even with sparsely sampled SDFs. Our approach provides a more nuanced understanding of SDFs and offers significant improvements in surface reconstruction.

Standard liquid simulators apply operator splitting to independently solve for pressure and viscous stresses, a decoupling that induces incorrect free surface boundary conditions. Such methods are unable to simulate fluid phenomena reliant on the balance of pressure and viscous stresses, such as the liquid rope coil instability exhibited by honey. By contrast, unsteady Stokes solvers retain coupling between pressure and viscosity, thus resolving these phenomena, albeit using a much larger and thus more computationally expensive linear system compared to the decoupled approach. To accelerate solving the unsteady Stokes problem, we propose a reduced fluid model wherein interior regions are represented with incompressible polynomial vector fields. Sets of standard grid cells are consolidated into super-cells, each of which are modelled using a quadratic field of 26 degrees of freedom. We demonstrate that the reduced field must necessarily be at least quadratic, with the affine model being unable to correctly capture viscous forces. We reproduce the liquid rope coiling instability, as well as other simulated examples, to show that our reduced model is able to reproduce the same fluid phenomena at a smaller computational cost. Futhermore, we performed a crowdsourced user survey to verify that our method produces imperceptible differences compared to the full unsteady Stokes method.

We introduce the walk-on-boundary (WoB) method for solving boundary value problems to computer graphics. WoB is a grid-free Monte Carlo solver for certain classes of second order partial differential equations. A similar Monte Carlo solver, the walk-on-spheres (WoS) method, has been recently popularized in computer graphics due to its advantages over traditional spatial discretization-based alternatives. We show that WoB's intrinsic properties yield further advantages beyond those of WoS. Unlike WoS, WoB naturally supports various boundary conditions (Dirichlet, Neumann, Robin, and mixed) for both interior and exterior domains. WoB builds upon boundary integral formulations, and it is mathematically more similar to light transport simulation in rendering than the random walk formulation of WoS. This similarity between WoB and rendering allows us to implement WoB on top of Monte Carlo ray tracing, and to incorporate advanced rendering techniques (e.g., bidirectional estimators with multiple importance sampling, the virtual point lights method, and Markov chain Monte Carlo) into WoB. WoB does not suffer from the intrinsic bias of WoS near the boundary and can estimate solutions precisely on the boundary. Our numerical results highlight the advantages of WoB over WoS as an attractive alternative to solve boundary value problems based on Monte Carlo.

Abstract We propose a novel surface‐only method for simulating dynamic deformables without the need for volumetric meshing or volumetric integral evaluations. While based upon a boundary element method (BEM)for linear elastodynamics, our method goes beyond simple adoption of BEM by addressing several of its key limitations. We alleviate large displacement artifacts due to linear elasticity by extending BEM with a moving reference frame and surface‐only fictitious forces, so that it only needs to handle deformations. To reduce memory and computational costs, we present a simple and practical method to compress the series of dense matrices required to simulate propagation of elastic waves over time. Furthermore, we explore a constraint enforcement mechanism and demonstrate the applicability of our method to general computer animation problems, such as frictional contact.

We propose to augment standard grid-based fluid solvers with pointwise divergence-free velocity interpolation, thereby ensuring exact incompressibility down to the sub-cell level. Our method takes as input a discretely divergence-free velocity field generated by a staggered grid pressure projection, and first recovers a corresponding discrete vector potential. Instead of solving a costly vector Poisson problem for the potential, we develop a fast parallel sweeping strategy to find a candidate potential and apply a gauge transformation to enforce the Coulomb gauge condition and thereby make it numerically smooth. Interpolating this discrete potential generates a point-wise vector potential whose analytical curl is a pointwise incompressible velocity field. Our method further supports irregular solid geometry through the use of level set-based cut-cells and a novel Curl-Noise-inspired potential ramping procedure that simultaneously offers strictly non-penetrating velocities and incompressibility. Experimental comparisons demonstrate that the vector potential reconstruction procedure at the heart of our approach is consistently faster than prior such reconstruction schemes, especially those that solve vector Poisson problems. Moreover, in exchange for its modest extra cost, our overall Curl-Flow framework produces significantly improved particle trajectories that closely respect irregular obstacles, do not suffer from spurious sources or sinks, and yield superior particle distributions over time.

We present a novel Monte Carlo-based fluid simulation approach capable of pointwise and stochastic estimation of fluid motion. Drawing on the Feynman-Kac representation of the vorticity transport equation, we propose a recursive Monte Carlo estimator of the Biot-Savart law and extend it with a stream function formulation that allows us to treat free-slip boundary conditions using a Walk-on-Spheres algorithm. Inspired by the Monte Carlo literature in rendering, we design and compare variance reduction schemes suited to a fluid simulation context for the first time, show its applicability to complex boundary settings, and detail a simple and practical implementation with temporal grid caching. We validate the correctness of our approach via quantitative and qualitative evaluations - across a range of settings and domain geometries - and thoroughly explore its parameters' design space. Finally, we provide an in-depth discussion of several axes of future work building on this new numerical simulation modality.

Simultaneous coupling of diverse physical systems poses significant computational challenges in terms of speed, quality, and stability. Rather than treating all components with a single discretization methodology (e.g., smoothed particles, material point method, Eulerian grid, etc.) that is ill-suited to some components, our solver, ElastoMonolith , addresses three-way interactions among standard particle-in-cell-based viscous and inviscid fluids, Lagrangian mesh-based deformable bodies, and rigid bodies. While prior methods often treat some terms explicitly or in a decoupled fashion for efficiency, often at the cost of robustness or stability, we demonstrate the effectiveness of a strong coupling approach that expresses all of the relevant physics within one consistent and unified optimization problem, including fluid pressure and viscosity, elasticity of the deformables, frictional solid-solid contact, and solid-fluid interface conditions. We further develop a numerical solver to tackle this difficult optimization problem, incorporating projected Newton, an active set method, and a transformation of the inner linear system matrix to ensure symmetric positive definiteness. Our experimental evaluations show that our framework can achieve high quality coupling results that avoid artifacts such as volume loss, instability, sticky contacts, and spurious interpenetrations.

The comprehensive visual modeling of fluid motion has historically been a challenging task, due in no small part to the difficulties inherent in geometries that are non-manifold, open, or thin. Modern geometric cut-cell mesh generators have been shown to produce, both robustly and quickly, workable volumetric elements in the presence of these problematic geometries, and the resulting volumetric representation would seem to offer an ideal infrastructure with which to perform fluid simulations. However, cut-cell mesh elements are general polyhedra that often contain holes and are non-convex; it is therefore difficult to construct the explicit function spaces required to employ standard functional discretizations, such as the Finite Element Method. The Virtual Element Method (VEM) has recently emerged as a functional discretization that successfully operates with complex polyhedral elements through a weak formulation of its function spaces. We present a novel cut-cell fluid simulation framework that exactly represents boundary geometry during the simulation. Our approach enables, for the first time, detailed fluid simulation with "in-the-wild" obstacles, including ones that contain non-manifold parts, self-intersections, and extremely thin features. Our key technical contribution is the generalization of the Particle-In-Cell fluid simulation methodology to arbitrary polyhedra using VEM. Coupled with a robust cut-cell generation scheme, this produces a fluid simulation algorithm that can operate on previously infeasible geometries without requiring any additional mesh modification or repair.

We propose FrictionalMonolith , a monolithic pressure-friction-contact solver for more accurately, robustly, and efficiently simulating two-way interactions of rigid bodies with continuum granular materials or inviscid liquids. By carefully formulating the components of such systems within a single unified minimization problem, our solver can simultaneously handle unilateral incompressibility and implicit integration of friction for the interior of the continuum, frictional contact resolution among the rigid bodies, and mutual force exchanges between the continuum and rigid bodies. Our monolithic approach eliminates various problematic artifacts in existing weakly coupled approaches, including loss of volume in the continuum material, artificial drift and slip of the continuum at solid boundaries, interpenetrations of rigid bodies, and simulation instabilities. To efficiently handle this challenging monolithic minimization problem, we present a customized solver for the resulting quadratically constrained quadratic program that combines elements of staggered projections, augmented Lagrangian methods, inexact projected Newton, and active-set methods. We demonstrate the critical importance of a unified treatment and the effectiveness of our proposed solver in a range of practical scenarios.

Abstract In this paper, we propose and evaluate fast, scalable approaches for solving the linear complementarity problems (LCP) arising from the fluid pressure equations with separating solid boundary conditions. Specifically, we present a policy iteration method, a penalty method, and a modified multigrid method, and demonstrate that each is able to properly handle the desired boundary conditions. Moreover, we compare our proposed methods against existing approaches and show that our solvers are more efficient and exhibit better scaling behavior; that is, the number of iterations required for convergence is essentially independent of grid resolution, and thus they are faster at larger grid resolutions. For example, on a 2563 grid our multigrid method was 30 times faster than the prior multigrid method in the literature.

We propose a new adaptive liquid simulation framework that achieves highly detailed behavior with reduced implementation complexity. Prior work has shown that spatially adaptive grids are efficient for simulating large-scale liquid scenarios, but in order to enable adaptivity along the liquid surface these methods require either expensive boundary-conforming (re-)meshing or elaborate treatments for second order accurate interface conditions. This complexity greatly increases the difficulty of implementation and maintainability, potentially making it infeasible for practitioners. We therefore present new algorithms for adaptive simulation that are comparatively easy to implement yet efficiently yield high quality results. First, we develop a novel staggered octree Poisson discretization for free surfaces that is second order in pressure and gives smooth surface motions even across octree T-junctions, without a power/Voronoi diagram construction. We augment this discretization with an adaptivity-compatible surface tension force that likewise supports T-junctions. Second, we propose a moving least squares strategy for level set and velocity interpolation that requires minimal knowledge of the local tree structure while blending near-seamlessly with standard trilinear interpolation in uniform regions. Finally, to maximally exploit the flexibility of our new surface-adaptive solver, we propose several novel extensions to sizing function design that enhance its effectiveness and flexibility. We perform a range of rigorous numerical experiments to evaluate the reliability and limitations of our method, as well as demonstrating it on several complex high-resolution liquid animation scenarios.

We propose to enhance the capability of standard free-surface flow simulators with efficient support for immersed bubbles through two new models: constraint-based bubbles and affine fluid regions. Unlike its predecessors, our constraint-based model entirely dispenses with the need for advection or projection inside zero-density bubbles, with extremely modest additional computational overhead that is proportional to the surface area of all bubbles. This surface-only approach is easy to implement, realistically captures many familiar bubble behaviors, and even allows two or more distinct liquid bodies to correctly interact across completely unsimulated air. We augment this model with a per-bubble volume-tracking and correction framework to minimize the cumulative effects of gradual volume drift. To support bubbles with non-zero densities, we propose a novel reduced model for an irregular fluid region with a single pointwise incompressible affine vector field. This model requires only 11 interior velocity degrees of freedom per affine fluid region in 3D, and correctly reproduces buoyant, stationary, and sinking behaviors of a secondary fluid phase with non-zero density immersed in water. Since the pressure projection step in both the above schemes is a slightly modified Poisson-style system, we propose novel Multigrid-based preconditioners for Conjugate Gradients for fast numerical solutions of our new discretizations. Furthermore, we observe that by enforcing an incompressible affine vector field over a coalesced set of grid cells, our reduced model is effectively an irregular coarse super-cell. This offers a convenient and flexible adaptive coarsening strategy that integrates readily with the standard staggered grid approach for fluid simulation, yet supports coarsened regions that are arbitrary voxelized shapes, and provides an analytically divergence-free interior. We demonstrate its effectiveness with a new adaptive liquid simulator whose interior regions are coarsened into a mix of tiles with regular and irregular shapes.

We propose Monolith , a monolithic pressure-viscosity-contact solver for more accurately, robustly, and efficiently simulating non-trivial two-way interactions of rigid bodies with inviscid, viscous, or non-Newtonian liquids. Our solver simultaneously handles incompressibility and (optionally) implicit viscosity integration for liquids, contact resolution for rigid bodies, and mutual interactions between liquids and rigid bodies by carefully formulating these as a single unified minimization problem. This monolithic approach reduces or eliminates an array of problematic artifacts, including liquid volume loss, solid interpenetrations, simulation instabilities, artificial "melting" of viscous liquid, and incorrect slip at liquid-solid interfaces. In the absence of solid-solid friction, our minimization problem is a Quadratic Program (QP) with a symmetric positive definite (SPD) matrix and can be treated with a single Linear Complementarity Problem (LCP) solve. When friction is present, we decouple the unified minimization problem into two subproblems so that it can be effectively handled via staggered projections with alternating LCP solves. We also propose a complementary approach for non-Newtonian fluids which can be seamlessly integrated and addressed during the staggered projections. We demonstrate the critical importance of a contact-aware , unified treatment of fluid-solid coupling and the effectiveness of our proposed Monolith solver in a wide range of practical scenarios.

We propose a framework for simulating the complex dynamics of strands interacting with compressible, shear-dependent liquids, such as oil paint, mud, cream, melted chocolate, and pasta sauce. Our framework contains three main components: the strands modeled as discrete rods, the bulk liquid represented as a continuum (material point method), and a reduced-dimensional flow of liquid on the surface of the strands with detailed elastoviscoplastic behavior. These three components are tightly coupled together. To enable discrete strands interacting with continuum-based liquid, we develop models that account for the volume change of the liquid as it passes through strands and the momentum exchange between the strands and the liquid. We also develop an extended constraint-based collision handling method that supports cohesion between strands. Furthermore, we present a principled method to preserve the total momentum of a strand and its surface flow, as well as an analytic plastic flow approach for Herschel-Bulkley fluid that enables stable semi-implicit integration at larger time steps. We explore a series of challenging scenarios, involving splashing, shaking, and agitating the liquid which causes the strands to stick together and become entangled.

While pressure forces are often the bottleneck in (near-)inviscid fluid simulations, viscosity can impose orders of magnitude greater computational costs at lower Reynolds numbers. We propose an implicit octree finite difference discretization that significantly accelerates the solution of the free surface viscosity equations using adaptive staggered grids, while supporting viscous buckling and rotation effects, variable viscosity, and interaction with scripted moving solids. In experimental comparisons against regular grids, our method reduced the number of active velocity degrees of freedom by as much as a factor of 7.7 and reduced linear system solution times by factors between 3.8 and 9.4. We achieve this by developing a novel adaptive variational finite difference methodology for octrees and applying it to the optimization form of the viscosity problem. This yields a linear system that is symmetric positive definite by construction, unlike naive finite difference/volume methods, and much sparser than a hypothetical finite element alternative. Grid refinement studies show spatial convergence at first order in L ∞ and second order in L 1 , while the significantly smaller size of the octree linear systems allows for the solution of viscous forces at higher effective resolutions than with regular grids. We demonstrate the practical benefits of our adaptive scheme by replacing the regular grid viscosity step of a commercial liquid simulator (Houdini) to yield large speed-ups, and by incorporating it into an existing inviscid octree simulator to add support for viscous flows. Animations of viscous liquids pouring, bending, stirring, buckling, and melting illustrate that our octree method offers significant computational gains and excellent visual consistency with its regular grid counterpart.

Although geometry arising "in the wild" most often comes in the form of a surface representation, a plethora of geometrical and physical applications require the construction of volumetric embeddings either of the geometry itself or the domain surrounding it. Cartesian cut-cell-based mesh generation provides an attractive solution in which volumetric elements are constructed from the intersection of the input surface geometry with a uniform or adaptive hexahedral grid. This choice, especially common in computational fluid dynamics, has the potential to efficiently generate accurate, surface-conforming cells; unfortunately, current solutions are often slow, fragile, or cannot handle many common topological situations. We therefore propose a novel, robust cut-cell construction technique for triangle surface meshes that explicitly computes the precise geometry of the intersection cells, even on meshes that are open or non-manifold. Its fundamental geometric primitive is the intersection of an arbitrary segment with an axis-aligned plane. Beginning from the set of intersection points between triangle mesh edges and grid planes, our bottom-up approach robustly determines cut-edges, cut-faces, and finally cut-cells, in a manner designed to guarantee topological correctness. We demonstrate its effectiveness and speed on a wide range of input meshes and grid resolutions, and make the code available as open source.

We introduce a new advection scheme for fluid animation. Our main contribution is the use of long-term temporal changes in pressure to extend the commonly used semi-Lagrangian scheme further back along the time axis. Our algorithm starts by tracing sample points along a trajectory following the velocity field backwards in time for many steps. During this backtracing process, the pressure gradient along the path is integrated to correct the velocity of the current time step. We show that our method effectively suppresses numerical diffusion, retains small-scale vorticity, and provides better long-term kinetic energy preservation.

We propose a method for simulating the complex dynamics of partially and fully saturated woven and knit fabrics interacting with liquid, including the effects of buoyancy, nonlinear drag, pore (capillary) pressure, dripping, and convection-diffusion. Our model evolves the velocity fields of both the liquid and solid relying on mixture theory, as well as tracking a scalar saturation variable that affects the pore pressure forces in the fluid. We consider the porous microstructure implied by the fibers composing individual threads, and use it to derive homogenized drag and pore pressure models that faithfully reflect the anisotropy of fabrics. In addition to the bulk liquid and fabric motion, we derive a quasi-static flow model that accounts for liquid spreading within the fabric itself. Our implementation significantly extends standard numerical cloth and fluid models to support the diverse behaviors of wet fabric, and includes a numerical method tailored to cope with the challenging nonlinearities of the problem. We explore a range of fabric-water interactions to validate our model, including challenging animation scenarios involving splashing, wringing, and collisions with obstacles, along with qualitative comparisons against simple physical experiments.

We present an efficient and scalable octree-inspired fluid simulation framework with the flexibility to leverage adaptivity in any part of the computational domain, even when resolution transitions reach the free surface. Our methodology ensures symmetry, definiteness and second order accuracy of the discrete Poisson operator, and eliminates numerical and visual artifacts of prior octree schemes. This is achieved by adapting the operators acting on the octree's simulation variables to reflect the structure and connectivity of a power diagram , which recovers primal-dual mesh orthogonality and eliminates problematic T-junction configurations. We show how such operators can be efficiently implemented using a pyramid of sparsely populated uniform grids, enhancing the regularity of operations and facilitating parallelization. A novel scheme is proposed for encoding the topology of the power diagram in the neighborhood of each octree cell, allowing us to locally reconstruct it on the fly via a lookup table, rather than resorting to costly explicit meshing. The pressure Poisson equation is solved via a highly efficient, matrix-free multigrid preconditioner for Conjugate Gradient, adapted to the power diagram discretization. We use another sparsely populated uniform grid for high resolution interface tracking with a narrow band level set representation. Using the recently introduced SPGrid data structure, sparse uniform grids in both the power diagram discretization and our narrow band level set can be compactly stored and efficiently updated via streaming operations. Additionally, we present enhancements to adaptive level set advection, velocity extrapolation, and the fast marching method for redistancing. Our overall framework gracefully accommodates the task of dynamically adapting the octree topology during simulation. We demonstrate end-to-end simulations of complex adaptive flows in irregularly shaped domains, with tens of millions of degrees of freedom.

The diverse interactions between hair and liquid are complex and span multiple length scales, yet are central to the appearance of humans and animals in many situations. We therefore propose a novel multi-component simulation framework that treats many of the key physical mechanisms governing the dynamics of wet hair. The foundations of our approach are a discrete rod model for hair and a particle-in-cell model for fluids. To treat the thin layer of liquid that clings to the hair, we augment each hair strand with a height field representation. Our contribution is to develop the necessary physical and numerical models to evolve this new system and the interactions among its components. We develop a new reduced-dimensional liquid model to solve the motion of the liquid along the length of each hair, while accounting for its moving reference frame and influence on the hair dynamics. We derive a faithful model for surface tension-induced cohesion effects between adjacent hairs, based on the geometry of the liquid bridges that connect them. We adopt an empirically-validated drag model to treat the effects of coarse-scale interactions between hair and surrounding fluid, and propose new volume-conserving dripping and absorption strategies to transfer liquid between the reduced and particle-in-cell liquid representations. The synthesis of these techniques yields an effective wet hair simulator, which we use to animate hair flipping, an animal shaking itself dry, a spinning car wash roller brush dunked in liquid, and intricate hair coalescence effects, among several additional scenarios.

We propose a novel unsteady Stokes solver for coupled viscous and pressure forces in grid-based liquid animation which yields greater accuracy and visual realism than previously achieved. Modern fluid simulators treat viscosity and pressure in separate solver stages, which reduces accuracy and yields incorrect free surface behavior. Our proposed implicit variational formulation of the Stokes problem leads to a symmetric positive definite linear system that gives properly coupled forces, provides unconditional stability, and treats difficult boundary conditions naturally through simple volume weights. Surface tension and moving solid boundaries are also easily incorporated. Qualitatively, we show that our method recovers the characteristic rope coiling instability of viscous liquids and preserves fine surface details, while previous grid-based schemes do not. Quantitatively, we demonstrate that our method is convergent through grid refinement studies on analytical problems in two dimensions. We conclude by offering practical guidelines for choosing an appropriate viscous solver, based on the scenario to be animated and the computational costs of different methods.

Abstract Combining high‐resolution level set surface tracking with lower resolution physics is an inexpensive method for achieving highly detailed liquid animations. Unfortunately, the inherent resolution mismatch introduces several types of disturbing visual artifacts. We identify the primary sources of these artifacts and present simple, efficient, and practical solutions to address them. First, we propose an unconditionally stable filtering method that selectively removes sub‐grid surface artifacts not seen by the fluid physics, while preserving fine detail in dynamic splashing regions. It provides comparable results to recent error‐correction techniques at lower cost, without substepping, and with better scaling behavior. Second, we show how a modified narrow‐band scheme can ensure accurate free surface boundary conditions in the presence of large resolution mismatches. Our scheme preserves the efficiency of the narrow‐band methodology, while eliminating objectionable stairstep artifacts observed in prior work. Third, we demonstrate that the use of linear interpolation of velocity during advection of the high‐resolution level set surface is responsible for visible grid‐aligned kinks; we therefore advocate higher‐order velocity interpolation, and show that it dramatically reduces this artifact. While these three contributions are orthogonal, our results demonstrate that taken together they efficiently address the dominant sources of visual artifacts arising with high‐resolution embedded liquid surfaces; the proposed approach offers improved visual quality, a straightforward implementation, and substantially greater scalability than competing methods.

Fluid animation methods based on Eulerian grids have long struggled to resolve flows involving narrow gaps and thin solid features. Past approaches have artificially inflated or voxelized boundaries, although this sacrifices the correct geometry and topology of the fluid domain and prevents flow through narrow regions. We present a boundary-respecting fluid simulator that overcomes these challenges. Our solution is to intersect the solid boundary geometry with the cells of a background regular grid to generate a topologically correct, boundary-conforming cut-cell mesh. We extend both pressure projection and velocity advection to support this enhanced grid structure. For pressure projection, we introduce a general graph-based scheme that properly preserves discrete incompressibility even in thin and topologically complex flow regions, while nevertheless yielding symmetric positive definite linear systems. For advection, we exploit polyhedral interpolation to improve the degree to which the flow conforms to irregular and possibly non-convex cell boundaries, and propose a modified PIC/FLIP advection scheme to eliminate the need to inaccurately reinitialize invalid cells that are swept over by moving boundaries. The method naturally extends the standard Eulerian fluid simulation framework, and while we focus on thin boundaries, our contributions are beneficial for volumetric solids as well. Our results demonstrate successful one-way fluid-solid coupling in the presence of thin objects and narrow flow regions even on very coarse grids.

We propose a novel surface-only technique for simulating incompressible, inviscid and uniform-density liquids with surface tension in three dimensions. The liquid surface is captured by a triangle mesh on which a Lagrangian velocity field is stored. Because advection of the velocity field may violate the incompressibility condition, we devise an orthogonal projection technique to remove the divergence while requiring the evaluation of only two boundary integrals. The forces of surface tension, gravity, and solid contact are all treated by a boundary element solve, allowing us to perform detailed simulations of a wide range of liquid phenomena, including waterbells, droplet and jet collisions, fluid chains, and crown splashes.

Simulating the delightful dynamics of soap films, bubbles, and foams has traditionally required the use of a fully three-dimensional many-phase Navier-Stokes solver, even though their visual appearance is completely dominated by the thin liquid surface. We depart from earlier work on soap bubbles and foams by noting that their dynamics are naturally described by a Lagrangian vortex sheet model in whichcirculationis the primary variable. This leads us to derive a novel circulation-preserving surface-only discretization of foam dynamics driven by surface tension on a non-manifold triangle mesh. We represent the surface using a mesh-based multimaterial surface tracker which supports complex bubble topology changes, and evolve the surface according to the ambient air flow induced by a scalar circulation field stored on the mesh. Surface tension forces give rise to a simple update rule for circulation, even at non-manifold Plateau borders, based on a discrete measure of signed scalar mean curvature. We further incorporate vertex constraints to enable the interaction of soap films with wires. The result is a method that is at once simple, robust, and efficient, yet able to capture an array of soap films behaviors including foam rearrangement, catenoid collapse, blowing bubbles, and double bubbles being pulled apart.

We consider the simulation of dense foams composed of microscopic bubbles, such as shaving cream and whipped cream. We represent foam not as a collection of discrete bubbles, but instead as a continuum. We employ thematerial point method(MPM) to discretize a hyperelastic constitutive relation augmented with the Herschel-Bulkley model of non-Newtonian viscoplastic flow, which is known to closely approximate foam behavior. Since large shearing flows in foam can produce poor distributions of material points, a typical MPM implementation can produce non-physical internal holes in the continuum. To address these artifacts, we introduce a particle resampling method for MPM. In addition, we introduce an explicit tearing model to prevent regions from shearing into artificially thin, honey-like threads. We evaluate our method's efficacy by simulating a number of dense foams, and we validate our method by comparing to real-world footage of foam.

We present a triangle mesh-based technique for tracking the evolution of three-dimensionalmultimaterial interfacesundergoing complex deformations. It is the first non-manifold triangle mesh tracking method to simultaneously maintain intersection-free meshes and support the proposed broad set of multimaterial remeshing and topological operations. We represent the interface as a non-manifold triangle mesh with material labels assigned to each half-face to distinguish volumetric regions. Starting from proposed application-dependent vertex velocities, we deform the mesh, seeking a non-intersecting, watertight solution. This goal necessitates development of various collision-safe, label-aware non-manifold mesh operations: multimaterial mesh improvement; T1 and T2 processes, topological transitions arising in foam dynamics and multiphase flows; and multimaterial merging, in which a new interface is created between colliding materials. We demonstrate the robustness and effectiveness of our approach on a range of scenarios including geometric flows and multiphase fluid animation.

We present the first reduced-dimensional technique to simulate the dynamics of thin sheets of viscous incompressible liquid in three dimensions. Beginning from a discrete Lagrangian model for elastic thin shells, we apply the Stokes-Rayleigh analogy to derive a simple yet consistent model for viscous forces. We incorporate nonlinear surface tension forces with a formulation based on minimizing discrete surface area, and preserve the quality of triangular mesh elements through local remeshing operations. Simultaneously, we track and evolve the thickness of each triangle to exactly conserve liquid volume. This approach enables the simulation of extremely thin sheets of viscous liquids, which are difficult to animate with existing volumetric approaches. We demonstrate our method with examples of several characteristic viscous sheet behaviors, including stretching, buckling, sagging, and wrinkling.

Abstract When simulating fluids, tetrahedral methods provide flexibility and ease of adaptivity that Cartesian grids find difficult to match. However, this approach has so far been limited by two conflicting requirements. First, accurate simulation requires quality Delaunay meshes and the use of circumcentric pressures. Second, meshes must align with potentially complex moving surfaces and boundaries, necessitating continuous remeshing. Unfortunately, sacrificing mesh quality in favour of speed yields inaccurate velocities and simulation artifacts. We describe how to eliminate the boundary‐matching constraint by adapting recent embedded boundary techniques to tetrahedra, so that neither air nor solid boundaries need to align with mesh geometry. This enables the use of high quality, arbitrarily graded, non‐conforming Delaunay meshes, which are simpler and faster to generate. Temporal coherence can also be exploited by reusing meshes over adjacent timesteps to further reduce meshing costs. Lastly, our free surface boundary condition eliminates the spurious currents that previous methods exhibited for slow or static scenarios. We provide several examples demonstrating that our efficient tetrahedral embedded boundary method can substantially increase the flexibility and accuracy of adaptive Eulerian fluid simulation.

We introduce an Eulerian liquid simulation framework based on the Voronoi diagram of a potentially unorganized collection of pressure samples. Constructing the simulation mesh in this way allows us to place samples anywhere in the computational domain; we exploit this by choosing samples that accurately capture the geometry and topology of the liquid surface. When combined with high-resolution explicit surface tracking this allows us to simulate nearly arbitrarily thin features, while eliminating noise and other artifacts that arise when there is a resolution mismatch between the simulation and the surface---and allowing a precise inclusion of surface tension based directly on and at the same resolution as the surface mesh. In addition, we present a simplified Voronoi/Delaunay mesh velocity interpolation scheme, and a direct extension of embedded free surfaces and solid boundaries to Voronoi meshes.

Physical simulation has emerged as a compelling animation technique, yet current approaches to coupling simulations of fluids and solids with irregular boundary geometry are inefficient or cannot handle some relevant scenarios robustly. We propose a new variational approach which allows robust and accurate solution on relatively coarse Cartesian grids, allowing possibly orders of magnitude faster simulation. By rephrasing the classical pressure projection step as a kinetic energy minimization, broadly similar to modern approaches to rigid body contact, we permit a robust coupling between fluid and arbitrary solid simulations that always gives a well-posed symmetric positive semi-definite linear system. We provide several examples of efficient fluid-solid interaction and rigid body coupling with sub-grid cell flow. In addition, we extend the framework with a new boundary condition for free-surface flow, allowing fluid to separate naturally from solids.

This article introduces the Hierarchical Run-Length Encoded (H-RLE) Level Set data structure. This novel data structure combines the best features of the DT-Grid (of Nielsen and Museth [2004]) and the RLE Sparse Level Set (of Houston et al. [2004]) to provide both optimal efficiency and extreme versatility. In brief, the H-RLE level set employs an RLE in a dimensionally recursive fashion. The RLE scheme allows the compact storage of sequential nonnarrowband regions while the dimensionally recursive encoding along each axis efficiently compacts nonnarrowband planes and volumes. Consequently, this new structure can store and process level sets with effective voxel resolutions exceeding 5000 × 3000 × 3000 (45 billion voxels) on commodity PCs with only 1 GB of memory. This article, besides introducing the H-RLE level set data structure and its efficient core algorithms, also describes numerous applications that have benefited from our use of this structure: our unified implicit object representation, efficient and robust mesh to level set conversion, rapid ray tracing, level set metamorphosis, collision detection, and fully sparse fluid simulation (including RLE vector and matrix representations.) Our comparisons of the popular octree level set and Peng level set structures to the H-RLE level set indicate that the latter is superior in both narrowband sequential access speed and overall memory usage.