Simulating heterogeneous multibody systems with both deformable and stiff components remains a challenge for GPU solvers. While Newton-Krylov methods are popular for their matrix-free nature and good GPU compatibility, they often suffer from severe ill-conditioning and slow convergence when handling stiff contacts and material disparities in deformable multibody systems. In this paper, we present Heterogeneous Subspace Corrections (HSC), a novel Newton-CG variant to efficiently simulate such complex systems. HSC decouples the system into two Krylov iterations: a Newton-CG procedure for deformable bodies and a Neumann-based iteration for the affine subsystem. We introduce a GPU-based Adaptive Cross Approximation (ACA) algorithm to exploit the low-rank nature of the coupling matrices, which manages to reduce the overhead of sparse matrix-vector multiplications substantially on the GPU. For the affine subsystem, we propose a dedicated data structure for fast assembling contact Hessians in parallel. A variety of experimental results demonstrate that our framework consistently outperforms existing GPU simulators, achieving convergence rates comparable to direct solvers even in scenes with high-resolution models and extensive stiff contacts.

We present an efficient B-spline finite element method (FEM) for cloth simulation. While higher-order FEM has long promised higher accuracy, its adoption in cloth simulators has been limited by its larger computational costs while generating results with similar visual quality. Our contribution is a full algorithmic pipeline that makes cloth simulation using quadratic B-spline surfaces faster than standard linear FEM in practice while consistently improving accuracy and visual fidelity. Using quadratic B-spline basis functions, we obtain a globally C 1 -continuous displacement field that supports consistent discretization of both membrane and bending energies, effectively reducing locking artifacts and mesh dependence common to linear elements. To close the performance gap, we introduce a reduced integration scheme that separately optimizes quadrature rules for membrane and bending energies, an accelerated Hessian assembly procedure tailored to the spline structure, and an optimized linear solver based on partial factorization. Together, these optimizations make high-order, smooth cloth simulation competitive at scale, yielding an average 2× speedup over linear FEM in our tests. Extensive experiments demonstrate improved accuracy, wrinkle detail, and robustness, including contact-rich scenarios, relative to linear FEM and recent higher-order approaches. Our method enables realistic wrinkling dynamics across a wide range of material parameters and supports practical garment animation, providing a new promising spatial discretization for high-quality cloth simulation.

Self-intersections are widespread in surface meshes and invalidate downstream simulation, fabrication, and learning pipelines. Existing approaches typically treat self-intersections as local collision events, but embeddedness (i.e., lack of self-intersections) is a global geometric property that cannot be enforced through local reasoning alone. We introduce an energy-based framework that enforces surface embeddedness simultaneously at the shape and mesh levels, based on the insight that successful untangling requires accounting for both global shape-level interactions and local mesh-level interactions. A shape-level energy captures global entanglement independent of discretization, while a mesh-level penalty regularizes local discrete interactions. Together, these energies enable reliable removal of self-intersections without changing mesh connectivity and apply to a broad class of geometries, including surfaces with boundary, non-manifold configurations, immersion failures, and multi-object scenes. Compared to prior state-of-the-art methods, our approach resolves self-intersections across challenging datasets, enabling reliable downstream processing of surface meshes.

Incremental Potential Contact (IPC) has emerged as a robust and unifying formulation for contact-rich physical simulation by casting elasticity and collision handling as a single energy minimization problem. Achieving high performance, however, typically requires heavily specialized implementations that hard-code assumptions about energies, primitive types, and parameterizations, creating a major barrier to extensibility. Adding new energies or alternative parameterizations often requires re-deriving first and second-order derivatives, and implementing new assembly logic for the global Hessian and gradient. This challenge is further exacerbated by collision energies, where the same energy definition is often applied to mixed parameterizations, which can lead to a combinatorial explosion of parameterization-specific derivative and assembly cases. In this paper we introduce YASPS, a framework for physical simulation that resolves this limitation by making structural relationships explicit in a differentiable representation. YASPS introduces two relational operators, JOIN and UNION, which encode connectivity and heterogeneous parameterizations directly in the symbolic computation graph. Using symbolic differentiation over these operators, YASPS automatically derives local derivatives, and determines the induced sparsity and block structure of global gradients and Hessians from the same description while avoiding any code explosion induced by mixed-parameterizations. Targeting IPC workloads, YASPS compiles the resulting symbolic graphs into GPU kernels for local energy evaluation, derivative computation, and block-sparse matrix assembly, and solves the resulting Newton systems using a GPU-based iterative solver. This approach achieves performance competitive with state-of-the-art IPC implementations while enabling new energies and parameterizations to be added through localized symbolic definitions, without hand-written derivative or assembly code.

We introduce a novel approach to simulate the interaction between fluids and thin elastic solids without any penetration. Our approach is centered around an optimization system augmented with barriers, which aims to find a configuration that ensures the absence of penetration while enforcing incompressibility for the fluids and minimizing elastic potentials for the solids. Unlike previous methods that primarily focus on velocity coherence at the fluid-solid interfaces, we demonstrate the effectiveness and flexibility of explicitly resolving positional constraints, including both explicit representation of solid positions and the implicit representation of fluid level-set interface. To preserve the volume of the fluid, we propose a simple yet efficient approach that adjusts the associated level-set values. Additionally, we develop a distance metric capable of measuring the separation between an implicitly represented surface and a Lagrangian object of arbitrary codimension. By integrating the inertia, solid elastic potential, damping, barrier potential, and fluid incompressibility within a unified system, we are able to robustly simulate a wide range of processes involving fluid interactions with lower-dimensional objects such as shells and rods. These processes include topology changes, bouncing, splashing, sliding, rolling, floating, and more.

Incremental Potential Contact (IPC) is a widely used, robust, and accurate method for simulating complex frictional contact behaviors. However, achieving high efficiency remains a major challenge, particularly as material stiffness increases, which leads to slower Preconditioned Conjugate Gradient (PCG) convergence, even with state-of-the-art preconditioners. In this article, we propose a fully GPU-optimized IPC simulation framework capable of handling materials across a wide range of stiffnesses, delivering consistent high performance and scalability with up to 10× speedup over state-of-the-art GPU IPC methods. Our framework introduces three key innovations: (1) A novel connectivity-enhanced Multilevel Additive Schwarz (MAS) preconditioner on the GPU, designed to efficiently capture both stiff and soft elastodynamics and improve PCG convergence at a reduced preconditioning cost. (2) A C 2 -continuous cubic energy with an analytic eigensystem for inexact strain limiting, enabling more parallel-friendly simulations of stiff membranes, such as cloth, without membrane locking. (3) For extremely stiff behaviors where elastic waves are barely visible, we employ affine body dynamics (ABD) with a hash-based two-level reduction strategy for fast Hessian assembly and efficient affine-deformable coupling. We conduct extensive performance analyses and benchmark studies to compare our framework against state-of-the-art methods and alternative design choices. Our system consistently delivers the fastest performance across soft, stiff, and hybrid simulation scenarios, even in cases with high resolution, large deformations, and high-speed impacts.

The Material Point Method (MPM) has become a cornerstone of physics-based simulation, widely used in geomechanics and computer graphics for modeling phenomena such as granular flows, viscoelasticity, fracture mechanics, etc. Despite its versatility, the original MPM suffers from cell-crossing instabilities caused by discontinuities in particle-grid transfer kernels. Existing solutions mostly mitigate these issues by adopting smoother shape functions, but at the cost of increased numerical diffusion and computational overhead due to larger kernel support. In this paper, we propose a novel C 2 -continuous compact kernel for MPM that achieves a unique balance in terms of stability, accuracy, and computational efficiency. Our method integrates seamlessly with Affine Particle-In-Cell (APIC) and Moving Least Squares (MLS) MPM, while only doubling the number of grid nodes associated with each particle compared to linear kernels. At its core is an innovative dual-grid framework, which associates particles with grid nodes exclusively within the cells they occupy on two staggered grids, ensuring consistent and stable force computations. We demonstrate that our method can be conveniently implemented using a domain-specific language, Taichi, or based on open-source GPU MPM frameworks, achieving faster runtime and less numerical diffusion compared to quadratic B-spline MPM. Comprehensive validation through unit tests, comparative studies, and stress tests demonstrates the efficacy of our approach in conserving both linear and angular momentum, handling stiff materials, and scaling efficiently for large-scale simulations. Our results highlight the transformative potential of compact, high-order kernels in advancing MPM's capabilities for stable, accurate, and high-performance simulations.

We introduce VR-Doh, an open-source, hands-on 3D modeling system that enables intuitive creation and manipulation of elastoplastic objects in Virtual Reality (VR). By customizing the Material Point Method (MPM) for real-time simulation of hand-induced large deformations and enhancing 3D Gaussian Splatting for seamless rendering, VR-Doh provides an interactive and immersive 3D modeling experience. Users can naturally sculpt, deform, and edit objects through both contact- and gesture-based hand-object interactions. To achieve real-time performance, our system incorporates localized simulation techniques, particle-level collision handling, and the decoupling of physical and appearance representations, ensuring smooth and responsive interactions. VR-Doh supports both object creation and editing, enabling diverse modeling tasks such as designing food items, characters, and interlocking structures, all resulting in simulation-ready assets. User studies with both novice and experienced participants highlight the system's intuitive design, immersive feedback, and creative potential. Compared to existing geometric modeling tools, VR-Doh offers enhanced accessibility and natural interaction, making it a powerful tool for creative exploration in VR.

We propose a GPU-based iterative method for accelerated elastodynamic simulation with the log-barrier-based contact model. While Newton's method is a conventional choice for solving the interior-point system, the presence of ill-conditioned log barriers often necessitates a direct solution at each linearized substep and costs substantial storage and computational overhead. Moreover, constraint sets that vary in each iteration present additional challenges in algorithm convergence. Our method employs a novel barrier-augmented Lagrangian method to improve system conditioning and solver efficiency by adaptively updating an augmentation constraint sets. This enables the utilization of a scalable, inexact Newton-PCG solver with sparse GPU storage, eliminating the need for direct factorization. We further enhance PCG convergence speed with a domain-decomposed warm start strategy based on an eigenvalue spectrum approximated through our in-time assembly. Demonstrating significant scalability improvements, our method makes simulations previously impractical on 128 GB of CPU memory feasible with only 8 GB of GPU memory and orders-of-magnitude faster. Additionally, our method adeptly handles stiff problems, surpassing the capabilities of existing GPU-based interior-point methods. Our results, validated across various complex collision scenarios involving intricate geometries and large deformations, highlight the exceptional performance of our approach.

We introduce a variational formulation for simulating sticky interactions between elastoplastic solids. Our method brings a wider range of material behaviors into the reach of the Incremental Potential Contact (IPC) solver recently developed by (Li et al. 2020). Extending IPC requires several contributions. We first augment IPC with the classical Raous-Cangemi-Cocou (RCC) adhesion model. This allows us to robustly simulate the sticky interactions between arbitrary codimensional-0, 1, and 2 geometries. To enable user-friendly practical adoptions of our method, we further introduce a physically parametrized, easily controllable normal adhesion formulation based on the unsigned distance, which is fully compatible with IPC's barrier formulation. Furthermore, we propose a smoothly clamped tangential adhesion model that naturally models intricate behaviors including debonding. Lastly, we perform benchmark studies comparing our method with the classical models as well as real-world experimental results to demonstrate the efficacy of our method.

Learning the long-range interactions on large-scale mesh-based physical systems with flat Graph Neural Networks (GNNs) and stacking Message Passings (MPs) is challenging due to the scaling complexity w.r.t. the number of nodes and over-smoothing. Therefore, there has been growing interest in the community to introduce multi-scale structures to GNNs for physics simulation. However, current state-of-the-art methods are limited by their reliance on the labor-heavy drawing of coarser meshes or building coarser levels based on spatial proximity, which can introduce wrong edges across geometry boundaries. Inspired by the bipartite graph determination, we propose a novel pooling strategy, bi-stride to tackle the aforementioned limitations. Bi-stride pools nodes on every other frontier of the Breadth-First-Search (BFS), without the need for the manual drawing of coarser meshes and, avoid wrong edges introduced by spatial proximity. Additionally, it enables a reduced number of MP times on each level and the non-parametrized pooling and unpooling by interpolations, similar to convolutional Neural Networks (CNNs), which significantly reduces computational requirements. Experiments show that the proposed framework, BSMS-GNN, significantly outperforms existing methods in terms of both accuracy and computational efficiency in representative physics-based simulation scenarios.

We propose an efficient cloth simulation method that combines the merits of two drastically different numerical procedures, namely the subspace integration and parallelizable iterative relaxation. We show those two methods can be organically coupled within the framework of projective dynamics (PD), where both low- and high-frequency cloth motions are effectively and efficiently computed. Our method works seamlessly with the state-of-the-art contact handling algorithm, the incremental potential contact (IPC), to offer the non-penetration guarantee of the resulting animation. Our core ingredient centers around the utilization of subspace for the expedited convergence of Jacobi-PD. This involves solving the reduced global system and smartly employing its precomputed factorization. In addition, we incorporate a time-splitting strategy to handle the frictional self-contacts.

We propose a hybrid neural network and physics framework for reduced-order modeling of elastoplasticity and fracture. State-of-the-art scientific computing models like the Material Point Method (MPM) faithfully simulate large-deformation elastoplasticity and fracture mechanics. However, their long runtime and large memory consumption render them unsuitable for applications constrained by computation time and memory usage, e.g., virtual reality. To overcome these barriers, we propose a reduced-order framework. Our key innovation is training a low-dimensional manifold for the Kirchhoff stress field via an implicit neural representation. This low-dimensional neural stress field (NSF) enables efficient evaluations of stress values and, correspondingly, internal forces at arbitrary spatial locations. In addition, we also train neural deformation and affine fields to build low-dimensional manifolds for the deformation and affine momentum fields. These neural stress, deformation, and affine fields share the same low-dimensional latent space, which uniquely embeds the high-dimensional simulation state. After training, we run new simulations by evolving in this single latent space, which drastically reduces the computation time and memory consumption. Our general continuum-mechanics-based reduced-order framework is applicable to any phenomena governed by the elastodynamics equation. To showcase the versatility of our framework, we simulate a wide range of material behaviors, including elastica, sand, metal, non-Newtonian fluids, fracture, contact, and collision. We demonstrate dimension reduction by up to 100,000 × and time savings by up to 10 ×.

In this paper, we present a GPU algorithm for finite element hyperelastic simulation. We show that the interior-point method, known to be effective for robust collision resolution, can be coupled with non-Newton procedures and be massively sped up on the GPU. Newton's method has been widely chosen for the interior-point family, which fully solves a linear system at each step. After that, the active set associated with collision/contact constraints is updated. Mimicking this routine using a non-Newton optimization (like gradient descent or ADMM) unfortunately does not deliver expected accelerations. This is because the barrier functions employed in an interior-point method need to be updated at every iteration to strictly confine the search to the feasible region. The associated cost (e.g., per-iteration CCD) quickly overweights the benefit brought by the GPU, and a new parallelism modality is needed. Our algorithm is inspired by the domain decomposition method and designed to move interior-point-related computations to local domains as much as possible. We minimize the size of each domain (i.e., a stencil) by restricting it to a single element, so as to fully exploit the capacity of modern GPUs. The stencil-level results are integrated into a global update using a novel hybrid sweep scheme. Our algorithm is locally second-order offering better convergence. It enables simulation acceleration of up to two orders over its CPU counterpart. We demonstrate the scalability, robustness, efficiency, and quality of our algorithm in a variety of simulation scenarios with complex and detailed collision geometries.

In this article, we present a four-layer distributed simulation system and its adaptation to the Material Point Method (MPM). The system is built upon a performance portable C++ programming model targeting major High-Performance-Computing (HPC) platforms. A key ingredient of our system is a hierarchical block-tile-cell sparse grid data structure that is distributable to an arbitrary number of Message Passing Interface (MPI) ranks. We additionally propose strategies for efficient dynamic load balance optimization to maximize the efficiency of MPI tasks. Our simulation pipeline can easily switch among backend programming models, including OpenMP and CUDA, and can be effortlessly dispatched onto supercomputers and the cloud. Finally, we construct benchmark experiments and ablation studies on supercomputers and consumer workstations in a local network to evaluate the scalability and load balancing criteria. We demonstrate massively parallel, highly scalable, and gigascale resolution MPM simulations of up to 1.01 billion particles for less than 323.25 seconds per frame with 8 OpenSSH-connected workstations.

We present a novel hybrid Lagrangian/Eulerian method for simulating inelastic flows that generates high-quality particle distributions with adaptive volumes. At its core, our approach integrates an updated Lagrangian time discretization of continuum mechanics with the Power Particle-In-Cell geometric representation of deformable materials. As a result, we obtain material points described by optimized density kernels that precisely track the varying particle volumes both spatially and temporally. For efficient CFL-rate simulations, we also propose an implicit time integration for our system using a non-linear Gauss-Seidel solver inspired by X-PBD, viewing Eulerian nodal velocities as primal variables. We demonstrate the versatility of our method with simulations of mesoscale bubbles, sands, liquid, and foams.

We present a robust and efficient method for simulating Lagrangian solid-fluid coupling based on a new operator splitting strategy. We use variational formulations to approximate fluid properties and solid-fluid interactions, and introduce a unified two-way coupling formulation for SPH fluids and FEM solids using interior point barrier-based frictional contact. We split the resulting optimization problem into a fluid phase and a solid-coupling phase using a novel time-splitting approach with augmented contact proxies , and propose efficient custom linear solvers. Our technique accounts for fluids interaction with nonlinear hyperelastic objects of different geometries and codimensions, while maintaining an algorithmically guaranteed non-penetrating criterion. Comprehensive benchmarks and experiments demonstrate the efficacy of our method.

In this paper, we propose a neural network-based approach for learning to represent the behavior of plastic solid materials ranging from rubber and metal to sand and snow. Unlike elastic forces such as spring forces, these plastic forces do not result from the positional gradient of any potential energy, imposing great challenges on the stability and flexibility of their simulation. Our method effectively resolves this issue by learning a generalizable plastic energy whose derivative closely matches the analytical behavior of plastic forces. Our method, for the first time, enables the simulation of a wide range of arbitrary elasticity-plasticity combinations using time step-independent, unconditionally stable optimization-based time integrators. We demonstrate the efficacy of our method by learning and producing challenging 2D and 3D effects of metal, sand, and snow with complex dynamics.

This study presents a new method for modeling the interaction between compressible flow, shock waves, and deformable structures, emphasizing destructive dynamics. Extending advances in time-splitting compressible flow and the Material Point Methods (MPM), we develop a hybrid Eulerian and Lagrangian/Eulerian scheme for monolithic flow-structure interactions. We adopt the second-order WENO scheme to advance the continuity equation. To stably resolve deforming boundaries with sub-cell particles, we propose a blending treatment of reflective and passable boundary conditions inspired by the theory of porous media. The strongly coupled velocity-pressure system is discretized with a new mixed-order finite element formulation employing B-spline shape functions. Shock wave propagation, temperature/density-induced buoyancy effects, and topology changes in solids are unitedly captured.

We present a simulation framework for multibody dynamics via a universal variational integration. Our method naturally supports mixed rigid-deformables and mixed codimensional geometries, while providing guaranteed numerical convergence and accurate resolution of contact, friction, and a wide range of articulation constraints. We unify (1) the treatment of simulation degrees of freedom for rigid and soft bodies by formulating them both in terms of Lagrangian nodal displacements, (2) the handling of general linear equality joint constraints through an efficient change-of-variable strategy, (3) the enforcement of nonlinear articulation constraints based on novel distance potential energies, (4) the resolution of frictional contact between mixed dimensions and bodies with a variational Incremental Potential Contact formulation, and (5) the modeling of generalized restitution through semi-implicit Rayleigh damping. We conduct extensive unit tests and benchmark studies to demonstrate the efficacy of our method.

Simulating stiff materials in applications where deformations are either not significant or else can safely be ignored is a fundamental task across fields. Rigid body modeling has thus long remained a critical tool and is, by far, the most popular simulation strategy currently employed for modeling stiff solids. At the same time, rigid body methods continue to pose a number of well known challenges and trade-offs including intersections, instabilities, inaccuracies, and/or slow performances that grow with contact-problem complexity. In this paper we revisit the stiff body problem and present ABD, a simple and highly effective affine body dynamics framework, which significantly improves state-of-the-art for simulating stiff-body dynamics. We trace the challenges in rigid-body methods to the necessity of linearizing piecewise-rigid trajectories and subsequent constraints. ABD instead relaxes the unnecessary (and unrealistic) constraint that each body's motion be exactly rigid with a stiff orthogonality potential, while preserving the rigid body model's key feature of a small coordinate representation. In doing so ABD replaces piecewise linearization with piecewise linear trajectories. This, in turn, combines the best of both worlds: compact coordinates ensure small, sparse system solves, while piecewise-linear trajectories enable efficient and accurate constraint (contact and joint) evaluations. Beginning with this simple foundation, ABD preserves all guarantees of the underlying IPC model we build it upon, e.g., solution convergence, guaranteed non-intersection, and accurate frictional contact. Over a wide range and scale of simulation problems we demonstrate that ABD brings orders of magnitude performance gains (two- to three-orders on the CPU and an order more when utilizing the GPU, obtaining 10, 000× speedups) over prior IPC-based methods, while maintaining simulation quality and nonintersection of trajectories. At the same time ABD has comparable or faster timings when compared to state-of-the-art rigid body libraries optimized for performance without guarantees, and successfully and efficiently solves challenging simulation problems where both classes of prior rigid body simulation methods fail altogether.

We present a GPU algorithm for deformable simulation. Our method offers good computational efficiency and penetration-free guarantee at the same time, which are not common with existing techniques. The main idea is an algorithmic integration of projective dynamics (PD) and incremental potential contact (IPC). PD is a position-based simulation framework, favored for its robust convergence and convenient implementation. We show that PD can be employed to handle the variational optimization with the interior point method e.g., IPC. While conceptually straightforward, this requires a dedicated rework over the collision resolution and the iteration modality to avoid incorrect collision projection with improved numerical convergence. IPC exploits a barrier-based formulation, which yields an infinitely large penalty when the constraint is on the verge of being violated. This mechanism guarantees intersection-free trajectories of deformable bodies during the simulation, as long as they are apart at the rest configuration. On the downside, IPC brings a large amount of nonlinearity to the system, making PD slower to converge. To mitigate this issue, we propose a novel GPU algorithm named A-Jacobi for faster linear solve at the global step of PD. A-Jacobi is based on Jacobi iteration, but it better harvests the computation capacity on modern GPUs by lumping several Jacobi steps into a single iteration. In addition, we also re-design the CCD root finding procedure by using a new minimum-gradient Newton algorithm. Those saved time budgets allow more iterations to accommodate stiff IPC barriers so that the result is both realistic and collision-free. Putting together, our algorithm simulates complicated models of both solids and shells on the GPU at an interactive rate or even in real time.

In this paper, we propose Energetically Consistent Inelasticity (ECI), a new formulation for modeling and discretizing finite strain elastoplasticity/viscoelasticity in a way that is compatible with optimization-based time integrators. We provide an in-depth analysis for allowing plasticity to be implicitly integrated through an augmented strain energy density function. We develop ECI on the associative von-Mises J2 plasticity, the non-associative Drucker-Prager plasticity, and the finite strain viscoelasticity. We demonstrate the resulting scheme on both the Finite Element Method (FEM) and the Material Point Method (MPM). Combined with a custom Newton-type optimization integration scheme, our method enables simulating stiff and large-deformation inelastic dynamics of metal, sand, snow, and foam with larger time steps, improved stability, higher efficiency, and better accuracy than existing approaches.

This paper introduces a new weighting scheme for particle-grid transfers that generates hybrid Lagrangian/Eulerian fluid simulations with uniform particle distributions and precise volume control. At its core, our approach reformulates the construction of Power Particles [de Goes et al. 2015] by computing volume-constrained density kernels. We employ these optimized kernels as particle domains within the Generalized Interpolation Material Point method (GIMP) in order to incorporate Power Particles into the Particle-In-Cell framework, hence the name the Power Particle-In-Cell method. We address the construction of volume-constrained density kernels as a regularized optimal transportation problem and describe an iterative solver based on localized Gaussian convolutions that leads to a significant performance speedup compared to [de Goes et al. 2015]. We also present novel extensions for handling free surfaces and solid obstacles that bypass the need for cell clipping and ghost particles. We demonstrate the advantages of our transfer weights by improving hybrid schemes for fluid simulation such as the Fluid Implicit Particle (FLIP) method and the Affine Particle-In-Cell (APIC) method with volume preservation and robustness to varying particle-per-cell ratio, while retaining low numerical dissipation, conserving linear and angular momenta, and avoiding particle reseeding or post-process relaxations.

Unmanned aerial vehicles have been demonstrated successfully in a variety of tasks, including surveying and sampling tasks over large areas. These vehicles can take many forms. Quadrotors’ agility and ability to hover makes them well suited for navigating potentially tight spaces, while fixed wing aircraft are capable of efficient flight over long distances. Hybrid aerial vehicles (HAVs) attempt to achieve both of these benefits by exhibiting multiple modes; however, morphing HAVs typically require extra actuators which add mass, reducing both agility and efficiency. We propose a morphing HAV with folding wings that exhibits both a quadrotor and a fixed wing mode without requiring any extra actuation. This is achieved by leveraging the motion of a bistable mechanism at the center of the aircraft to drive folding of the wing using only the existing motors and the inertia of the system. We optimize both the bistable mechanism and the folding wing using a topology optimization approach. The resulting mechanisms were fabricated on a 3D printer and replaced the frame of an existing quadrotor. Our prototype successfully transitions between both modes and our experiments demonstrate that the behavior of the fabricated prototype is consistent with that of the simulation.

We extend recent advances in the numerical time-integration of contacting elastodynamics [Li et al. 2020] to build a new framework, called Injective Deformation Processing (IDP), for the robust solution of a wide range of mesh deformation problems requiring injectivity. IDP solves challenging 3D (and 2D) geometry processing and animation tasks on meshes, via artificial time integration, with guarantees of both non-inversion and non-overlap. To our knowledge IDP is the first framework for 3D deformation processing that can efficiently guarantee globally injective deformation without geometric locking. We demonstrate its application on a diverse set of problems and show its significant improvement over state-of-the-art for globally injective 3D deformation.

We introduce the first implicit time-stepping algorithm for rigid body dynamics, with contact and friction, that guarantees intersection-free configurations at every time step. Our algorithm explicitly models the curved trajectories traced by rigid bodies in both collision detection and response. For collision detection, we propose a conservative narrow phase collision detection algorithm for curved trajectories, which reduces the problem to a sequence of linear CCD queries with minimal separation. For time integration and contact response, we extend the recently proposed incremental potential contact framework to reduced coordinates and rigid body dynamics. We introduce a benchmark for rigid body simulation and show that our approach, while less efficient than alternatives, can robustly handle a wide array of complex scenes, which cannot be simulated with competing methods, without requiring per-scene parameter tuning.

We propose a framework of efficient nonlinear deformable simulation with both fast continuous collision detection and robust collision resolution. We name this new framework Medial IPC as it integrates the merits from medial elastics, for an efficient and versatile reduced simulation, as well as incremental potential contact, for a robust collision and contact resolution. We leverage medial axis transform to construct a kinematic subspace. Instead of resorting to projective dynamics, we use classic hyperelastics to embrace real-world nonlinear materials. A novel reduced continuous collision detection algorithm is presented based on the medial mesh. Thanks to unique geometric properties of medial axis and medial primitives, we derive closed-form formulations for identifying between-primitive collision within the reduced medial space. In the meantime, the implicit barrier energy that generates necessary repulsion forces for collision resolution is also formulated with the medial coordinate. In other words, Medial IPC exploits a universal reduced coordinate for simulation, continuous self-/collision detection, and IPC-based collision resolution. Continuous collision detection also allows more aggressive time stepping. In addition, we carefully implement our system with a heterogeneous CPU-GPU deployment such that massively parallelizable computations are carried out on the GPU while few sequential computations are on the CPU. Such implementation also frees us from generating training poses for selecting Cubature points and pre-computing their weights. We have tested our method on complicated deformable models and collision-rich simulation scenarios. Due to the reduced nature of our system, the computation is faster than fullspace IPC or other fullspace methods using continuous collision detection by at least one order. The simulation remains high-quality as the medial subspace captures intriguing and local deformations with sufficient realism.

We extend the incremental potential contact (IPC) model [Li et al. 2020a] for contacting elastodynamics to resolve systems composed of codimensional degrees-of-freedoms in arbitrary combination. This enables a unified, interpenetration-free, robust, and stable simulation framework that couples codimension-0,1,2, and 3 geometries seamlessly with frictional contact. Extending the IPC model to thin structures poses new challenges in computing strain, modeling thickness and determining collisions. To address these challenges we propose three corresponding contributions. First, we introduce a C 2 constitutive barrier model that directly enforces strain limiting as an energy potential while preserving rest state. This provides energetically-consistent strain limiting models (both isotropic and anisotropic) for cloth that enable strict satisfaction of strain-limit inequalities with direct coupling to both elastodynamics and contact via minimization of the incremental potential. Second, to capture the geometric thickness of codimensional domains we extend the IPC model to directly enforce distance offsets. Our treatment imposes a strict guarantee that mid-surfaces (respectively mid-lines) of shells (respectively rods) will not move closer than applied thickness values, even as these thicknesses become characteristically small. This enables us to account for thickness in the contact behavior of codimensional structures and so robustly capture challenging contacting geometries; a number of which, to our knowledge, have not been simulated before. Third, codimensional models, especially with modeled thickness, mandate strict accuracy requirements that pose a severe challenge to all existing continuous collision detection (CCD) methods. To address these limitations we develop a new, efficient, simple-to-implement additive CCD (ACCD) method that applies conservative advancement [Mirtich 1996; Zhang et al. 2006] to iteratively refine a lower bound for deforming primitives, converging to time of impact. In combination these contributions enable codimensional IPC (C-IPC). We perform extensive benchmark experiments to validate the efficacy of our method in capturing intricate behaviors of thin-structure contact and resulting bulk effects. In our experiments C-IPC obtains feasible, convergent, and so artifact-free solutions for all time steps, across all tested examples - producing robust simulations. We test C-IPC across extreme deformations, large time steps, and exceedingly close contact over all possible pairings of codimensional domains. Finally, with our strain-limit model, we confirm C-IPC guarantees non-intersection and strain-limit satisfaction for all reasonable (and well below - verified down to 0.1%) strain limits throughout all time steps.

We propose a novel scheme for simulating two-way coupled interactions between nonlinear elastic solids and incompressible fluids. The key ingredient of this approach is a ghost matrix operator-splitting scheme for strongly coupled nonlinear elastica and incompressible fluids through the weak form of their governing equations. This leads to a stable and efficient method handling large time steps under the CFL limit while using a single monolithic solve for the coupled pressure fields, even in the case with highly nonlinear elastic solids. The use of the Material Point Method (MPM) is essential in the designing of the scheme, it not only preserves discretization consistency with the hybrid Lagrangian-Eulerian fluid solver, but also works naturally with our novel interface quadrature (IQ) discretization for free-slip boundary conditions. While traditional MPM suffers from sticky numerical artifacts, our framework naturally supports discontinuous tangential velocities at the solid-fluid interface. Our IQ discretization results in an easy-to-implement, fully particle-based treatment of the interfacial boundary, avoiding the additional complexities associated with intermediate level set or explicit mesh representations. The efficacy of the proposed scheme is verified by various challenging simulations with fluid-elastica interactions.

Contacts weave through every aspect of our physical world, from daily household chores to acts of nature. Modeling and predictive computation of these phenomena for solid mechanics is important to every discipline concerned with the motion of mechanical systems, including engineering and animation. Nevertheless, efficiently time-stepping accurate and consistent simulations of real-world contacting elastica remains an outstanding computational challenge. To model the complex interaction of deforming solids in contact we propose Incremental Potential Contact (IPC) - a new model and algorithm for variationally solving implicitly time-stepped nonlinear elastodynamics. IPC maintains an intersection- and inversion-free trajectory regardless of material parameters, time step sizes, impact velocities, severity of deformation, or boundary conditions enforced. Constructed with a custom nonlinear solver, IPC enables efficient resolution of time-stepping problems with separate, user-exposed accuracy tolerances that allow independent specification of the physical accuracy of the dynamics and the geometric accuracy of surface-to-surface conformation. This enables users to decouple, as needed per application, desired accuracies for a simulation's dynamics and geometry. The resulting time stepper solves contact problems that are intersection-free (and thus robust), inversion-free, efficient (at speeds comparable to or faster than available methods that lack both convergence and feasibility), and accurate (solved to user-specified accuracies). To our knowledge this is the first implicit time-stepping method, across both the engineering and graphics literature that can consistently enforce these guarantees as we vary simulation parameters. In an extensive comparison of available simulation methods, research libraries and commercial codes we confirm that available engineering and computer graphics methods, while each succeeding admirably in custom-tuned regimes, often fail with instabilities, egregious constraint violations and/or inaccurate and implausible solutions, as we vary input materials, contact numbers and time step. We also exercise IPC across a wide range of existing and new benchmark tests and demonstrate its accurate solution over a broad sweep of reasonable time-step sizes and beyond (up to h = 2s) across challenging large-deformation, large-contact stress-test scenarios with meshes composed of up to 2.3M tetrahedra and processing up to 498K contacts per time step. For applications requiring high-accuracy we demonstrate tight convergence on all measures. While, for applications requiring lower accuracies, e.g. animation, we confirm IPC can ensure feasibility and plausibility even when specified tolerances are lowered for efficiency.

We propose Hierarchical Optimization Time Integration (HOT) for efficient implicit timestepping of the material point method (MPM) irrespective of simulated materials and conditions. HOT is an MPM-specialized hierarchical optimization algorithm that solves nonlinear timestep problems for large-scale MPM systems near the CFL limit. HOT provides convergent simulations out of the box across widely varying materials and computational resolutions without parameter tuning. As an implicit MPM timestepper accelerated by a custom-designed Galerkin multigrid wrapped in a quasi-Newton solver, HOT is both highly parallelizable and robustly convergent. As we show in our analysis, HOT maintains consistent and efficient performance even as we grow stiffness, increase deformation, and vary materials over a wide range of finite strain, elastodynamic, and plastic examples. Through careful benchmark ablation studies, we compare the effectiveness of HOT against seemingly plausible alternative combinations of MPM with standard multigrid and other Newton-Krylov models. We show how these alternative designs result in severe issues and poor performance. In contrast, HOT outperforms existing state-of-the-art, heavily optimized implicit MPM codes with an up to 10× performance speedup across a wide range of challenging benchmark test simulations.

Harnessing the power of modern multi-GPU architectures, we present a massively parallel simulation system based on the Material Point Method (MPM) for simulating physical behaviors of materials undergoing complex topological changes, self-collision, and large deformations. Our system makes three critical contributions. First, we introduce a new particle data structure that promotes coalesced memory access patterns on the GPU and eliminates the need for complex atomic operations on the memory hierarchy when writing particle data to the grid. Second, we propose a kernel fusion approach using a new Grid-to-Particles-to-Grid ( G2P2G ) scheme, which efficiently reduces GPU kernel launches, improves latency, and significantly reduces the amount of global memory needed to store particle data. Finally, we introduce optimized algorithmic designs that allow for efficient sparse grids in a shared memory context, enabling us to best utilize modern multi-GPU computational platforms for hybrid Lagrangian-Eulerian computational patterns. We demonstrate the effectiveness of our method with extensive benchmarks, evaluations, and dynamic simulations with elastoplasticity, granular media, and fluid dynamics. In comparisons against an open-source and heavily optimized CPU-based MPM codebase [Fang et al. 2019] on an elastic sphere colliding scene with particle counts ranging from 5 to 40 million, our GPU MPM achieves over 100x per-time-step speedup on a workstation with an Intel 8086K CPU and a single Quadro P6000 GPU, exposing exciting possibilities for future MPM simulations in computer graphics and computational science. Moreover, compared to the state-of-the-art GPU MPM method [Hu et al. 2019a], we not only achieve 2x acceleration on a single GPU but our kernel fusion strategy and Array-of-Structs-of-Array ( AoSoA ) data structure design also generalizes to multi-GPU systems. Our multi-GPU MPM exhibits near-perfect weak and strong scaling with 4 GPUs, enabling performant and large-scale simulations on a 1024 3 grid with close to 100 million particles with less than 4 minutes per frame on a single 4-GPU workstation and 134 million particles with less than 1 minute per frame on an 8-GPU workstation.

Dynamic fracture surrounds us in our day-to-day lives, but animating this phenomenon is notoriously difficult and only further complicated by anisotropic materials---those with underlying structures that dictate preferred fracture directions. Thus, we present AnisoMPM: a robust and general approach for animating the dynamic fracture of isotropic, transversely isotropic, and orthotropic materials. AnisoMPM has three core components: a technique for anisotropic damage evolution, methods for anisotropic elastic response, and a coupling approach. For anisotropic damage, we adopt a non-local continuum damage mechanics (CDM) geometric approach to crack modeling and augment this with structural tensors to encode material anisotropy. Furthermore, we discretize our damage evolution with explicit and implicit integration, giving a high degree of computational efficiency and flexibility. We also utilize a QR-decomposition based anisotropic constitutive model that is inversion safe, more efficient than SVD models, easy to implement, robust to extreme deformations, and that captures all aforementioned modes of anisotropy. Our elasto-damage coupling is enforced through an additive decomposition of our hyperelasticity into a tensile and compressive component in which damage is used to degrade the tensile contribution to allow for material separation. For extremely stiff fibered materials, we further introduce a novel Galerkin weak form discretization that enables embedded directional inextensibility. We present this as a hard-constrained grid velocity solve that poses an alternative to our anisotropic elasticity that is locking-free and can model very stiff materials.

Simulating viscoelastic polymers and polymeric fluids requires a robust and accurate capture of elasticity and viscosity. The computation is known to become very challenging under large deformations and high viscosity. Drawing inspirations from return mapping based elastoplasticity treatment for granular materials, we present a finite strain integration scheme for general viscoelastic solids under arbitrarily large deformation and non-equilibrated flow. Our scheme is based on a predictor-corrector exponential mapping scheme on the principal strains from the deformation gradient, which closely resembles the conventional treatment for elastoplasticity and allows straightforward implementation into any existing constitutive models. We develop a new Material Point Method that is fully implicit on both elasticity and inelasticity using augmented Lagrangian optimization with various preconditioning strategies for highly efficient time integration. Our method not only handles viscoelasticity but also supports existing elastoplastic models including Drucker-Prager and von-Mises in a unified manner. We demonstrate the efficacy of our framework on various examples showing intricate and characteristic inelastic dynamics with competitive performance.

Simulation methods are rapidly advancing the accuracy, consistency and controllability of elastodynamic modeling and animation. Critical to these advances, we require efficient time step solvers that reliably solve all implicit time integration problems for elastica. While available time step solvers succeed admirably in some regimes, they become impractically slow, inaccurate, unstable, or even divergent in others --- as we show here. Towards addressing these needs we present the Decomposed Optimization Time Integrator (DOT), a new domain-decomposed optimization method for solving the per time step, nonlinear problems of implicit numerical time integration. DOT is especially suitable for large time step simulations of deformable bodies with nonlinear materials and high-speed dynamics. It is efficient, automated, and robust at large, fixed-size time steps, thus ensuring stable, continued progress of high-quality simulation output. Across a broad range of extreme and mild deformation dynamics, using frame-rate size time steps with widely varying object shapes and mesh resolutions, we show that DOT always converges to user-set tolerances, generally well-exceeding and always close to the best wall-clock times across all previous nonlinear time step solvers, irrespective of the deformation applied.

We present two new approaches for animating dynamic fracture involving large elastoplastic deformation. In contrast to traditional mesh-based techniques, where sharp discontinuity is introduced to split the continuum at crack surfaces, our methods are based on Continuum Damage Mechanics (CDM) with a variational energy-based formulation for crack evolution. Our first approach formulates the resulting dynamic material damage evolution with a Ginzburg-Landau type phase-field equation and discretizes it with the Material Point Method (MPM), resulting in a coupled momentum/damage solver rooted in phase field fracture: PFF-MPM. Although our PFF-MPM approach achieves convincing fracture with or without plasticity, we also introduce a return mapping algorithm that can be analytically solved for a wide range of general non-associated plasticity models, achieving more than two times speedup over traditional iterative approaches. To demonstrate the efficacy of the algorithm, we also develop a Non-Associated Cam-Clay (NACC) plasticity model with a novel fracture-friendly hardening scheme. Our NACC plasticity paired with traditional MPM composes a second approach to dynamic fracture, as it produces a breadth of organic, brittle material fracture effects on its own. Though NACC and PFF can be combined, we focus on exploring their material effects separately. Both methods can be easily integrated into any existing MPM solver, enabling the simulation of various fracturing materials with extremely high visual fidelity while requiring little additional computational overhead.

Low-distortion mapping of three-dimensional surfaces to the plane is a critical problem in geometry processing. The intrinsic distortion introduced by these UV mappings is highly dependent on the choice of surface cuts that form seamlines which break mapping continuity. Parameterization applications typically require UV maps with an application-specific upper bound on distortion to avoid mapping artifacts; at the same time they seek to reduce cut lengths to minimize discontinuity artifacts. We propose OptCuts , an algorithm that jointly optimizes the parameterization and cutting of a three-dimensional mesh. OptCuts starts from an arbitrary initial embedding and a user-requested distortion bound. It requires no parameter setting and automatically seeks to minimize seam lengths subject to satisfying the distortion bound of the mapping computed using these seams. OptCuts alternates between topology and geometry update steps that consistently decrease distortion and seam length, producing a UV map with compact boundaries that strictly satisfies the distortion bound. OptCuts automatically produces high-quality, globally bijective UV maps without user intervention. While OptCuts can thus be a highly effective tool to create new mappings from scratch, we also show how it can be employed to improve pre-existing embeddings. Additionally, when semantic or other priors on seam placement are desired, OptCuts can be extended to respect these user preferences as constraints during optimization of the parameterization. We demonstrate the scalable performance of OptCuts on a wide range of challenging benchmark parameterization examples, as well as in comparisons with state-of-the-art UV methods and commercial tools.

While folds and pleats add interest to garments and cloth objects, incorporating them into an existing design manually or using existing software requires expertise and time. We present FoldSketch , a new system that supports simple and intuitive fold and pleat design. FoldSketch users specify the fold or pleat configuration they seek using a simple schematic sketching interface; the system then algorithmically generates both the fold-enhanced 3D garment geometry that conforms to user specifications, and the corresponding 2D patterns that reproduce this geometry within a simulation engine. While previous work aspired to compute the desired patterns for a given target 3D garment geometry, our main algorithmic challenge is that we do not have target geometry to start with. Real-life garment folds have complex profile shapes, and their exact geometry and location on a garment are intricately linked to a range of physical factors such as fabric properties and the garment's interaction with the wearer's body; it is therefore virtually impossible to predict the 3D shape of a fold-enhanced garment using purely geometric means. At the same time, using physical simulation to model folds requires appropriate 2D patterns and initial drape, neither of which can be easily provided by the user. We obtain both the 3D fold-enhanced garment and its corresponding patterns and initial drape via an alternating 2D-3D algorithm. We first expand the input patterns by allocating excess material for the expected fold formation; we then use these patterns to produce an estimated fold-enhanced drape geometry that balances designer expectations against physical reproducibility. We use the patterns and the estimated drape as input to a simulation generating an initial reproducible output. We improve the output's alignment with designer expectations by progressively refining the patterns and the estimated drape, converging to a final fully physically reproducible fold-enhanced garment. Our experiments confirm that FoldSketch reliably converges to a desired garment geometry and corresponding patterns and drape, and works well with different physical simulators. We demonstrate the versatility of our approach by showcasing a collection of garments augmented with diverse fold and pleat layouts specified via the FoldSketch interface, and further validate our approach via comparisons to alternative solutions and feedback from potential users.

Most fluid scenarios in graphics have a high Reynolds number, where viscosity is dominated by inertial effects, thus most solvers drop viscosity altogether: numerical damping from coarse grids is generally stronger than physical viscosity while resembling it in character. However, viscosity remains crucial near solid boundaries, in the boundary layer , to a large extent determining the look of the flow as a function of Reynolds number. Typical graphics simulations do not resolve boundary layer dynamics, so their look is determined mostly by numerical errors with the given grid size and time step, rather than physical parameters. We introduce two complementary techniques to capture boundary layer dynamics, bringing more physical control and predictability. We extend the FLIP particle-grid method with viscous particle strength exchange[Rivoalen and Huberson 2001] to better transfer momentum at solid boundaries, dubbed VFLIP. We also introduce Weakly Higher Resolution Regional Projection (WHIRP), a cheap and simple way to increase grid resolution where important by overlaying high resolution grids on the global coarse grid.

Recent studies have investigated a novel software paradigm that enables program decomposition and dynamic partitioning for cloud gaming system. Among the challenges for the proposed software architecture, to design and implement the decomposed software, which follows the platform's specification, is most critical in practice. In this work, we investigate a fundamental application in gaming scenario, the 3D skeletal game engine, which is often used to animate avatars. A graphics demo program is developed to demonstrate the feasibility of proposed system, while experiments have been conducted to show the efficiency of cognitive resource allocation. We also made contributions to seeking the appropriate design pattern in developing programs on our platform. In addition, we provide future visions on alternative implementations of some specific parts, which dedicates on improving the platform's compatibility.