Search

The eigenfunctions of the Laplace operator are essential in mathematical physics, engineering, and geometry processing. Typically, these are computed by discretizing the domain and performing eigendecomposition, tying the results to a specific mesh. However, this method is unsuitable for continuously-parameterized shapes. We propose a novel representation for eigenfunctions in continuously-parameterized shape spaces, where eigenfunctions are spatial fields with continuous dependence on shape parameters, defined by minimal Dirichlet energy, unit norm, and mutual orthogonality. We implement this with multilayer perceptrons trained as neural fields, mapping shape parameters and domain positions to eigenfunction values. A unique challenge is enforcing mutual orthogonality with respect to causality, where the causal ordering varies across the shape space. Our training method therefore requires three interwoven concepts: (1) learning $n$ eigenfunctions concurrently by minimizing Dirichlet energy with unit norm constraints; (2) filtering gradients during backpropagation to enforce causal orthogonality, preventing earlier eigenfunctions from being influenced by later ones; (3) dynamically sorting the causal ordering based on eigenvalues to track eigenvalue curve crossovers. We demonstrate our method on problems such as shape family analysis, predicting eigenfunctions for incomplete shapes, interactive shape manipulation, and computing higher-dimensional eigenfunctions, on all of which traditional methods fall short.

Traditional character animation specializes in characters with a rigidly articulated skeleton and a bipedal/quadripedal morphology. This assumption simplifies many aspects for designing physically based animations, like locomotion, but comes with the price of excluding characters of arbitrary deformable geometries. To remedy this, our framework makes use of a spatio-temporal actuation subspace built off of the natural vibration modes of the character geometry. The resulting actuation is coupled to a reduced fast soft body simulation, allowing us to formulate a locomotion optimization problem that is tractable for a wide variety of high resolution deformable characters., Comment: 15 pages, 14 figures

We propose a reduced space mixed finite element method (MFEM) built on a Skinning Eigenmode subspace and material-aware cubature scheme. Our solver is well-suited for simulating scenes with large material and geometric heterogeneities in real-time. This mammoth geometry is composed of 98,175 vertices and 531,565 tetrahedral elements and with a heterogenous composition of widely varying materials of muscles ($E= 5\times10^5$ Pa), joints ($E=1\times10^5$ Pa), and bone ($E=1\times10^{10}$ Pa). The resulting simulation runs at 120 frames per second (FPS)., Comment: 10 pages, 15 figures

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,000X and time savings by up to 10X.

Linear reduced-order modeling (ROM) simplifies complex simulations by approximating the behavior of a system using a simplified kinematic representation. Typically, ROM is trained on input simulations created with a specific spatial discretization, and then serves to accelerate simulations with the same discretization. This discretization-dependence is restrictive. Becoming independent of a specific discretization would provide flexibility to mix and match mesh resolutions, connectivity, and type (tetrahedral, hexahedral) in training data; to accelerate simulations with novel discretizations unseen during training; and to accelerate adaptive simulations that temporally or parametrically change the discretization. We present a flexible, discretization-independent approach to reduced-order modeling. Like traditional ROM, we represent the configuration as a linear combination of displacement fields. Unlike traditional ROM, our displacement fields are continuous maps from every point on the reference domain to a corresponding displacement vector; these maps are represented as implicit neural fields. With linear continuous ROM (LiCROM), our training set can include multiple geometries undergoing multiple loading conditions, independent of their discretization. This opens the door to novel applications of reduced order modeling. We can now accelerate simulations that modify the geometry at runtime, for instance via cutting, hole punching, and even swapping the entire mesh. We can also accelerate simulations of geometries unseen during training. We demonstrate one-shot generalization, training on a single geometry and subsequently simulating various unseen geometries.

We propose a reduced-space elasto-dynamic solver that is well suited for augmenting rigged character animations with secondary motion. At the core of our method is a novel deformation subspace based on Linear Blend Skinning that overcomes many of the shortcomings prior subspace methods face. Our skinning subspace is parameterized entirely by a set of scalar weights, which we can obtain through a small, material-aware and rig-sensitive generalized eigenvalue problem. The resulting subspace can easily capture rotational motion and guarantees that the resulting simulation is rotation equivariant. We further propose a simple local-global solver for linear co-rotational elasticity and propose a clustering method to aggregate per-tetrahedra non-linear energetic quantities. The result is a compact simulation that is fully decoupled from the complexity of the mesh., Comment: 20 pages, 24 figures

Implicit Neural Spatial Representation (INSR) has emerged as an effective representation of spatially-dependent vector fields. This work explores solving time-dependent PDEs with INSR. Classical PDE solvers introduce both temporal and spatial discretizations. Common spatial discretizations include meshes and meshless point clouds, where each degree-of-freedom corresponds to a location in space. While these explicit spatial correspondences are intuitive to model and understand, these representations are not necessarily optimal for accuracy, memory usage, or adaptivity. Keeping the classical temporal discretization unchanged (e.g., explicit/implicit Euler), we explore INSR as an alternative spatial discretization, where spatial information is implicitly stored in the neural network weights. The network weights then evolve over time via time integration. Our approach does not require any training data generated by existing solvers because our approach is the solver itself. We validate our approach on various PDEs with examples involving large elastic deformations, turbulent fluids, and multi-scale phenomena. While slower to compute than traditional representations, our approach exhibits higher accuracy and lower memory consumption. Whereas classical solvers can dynamically adapt their spatial representation only by resorting to complex remeshing algorithms, our INSR approach is intrinsically adaptive. By tapping into the rich literature of classic time integrators, e.g., operator-splitting schemes, our method enables challenging simulations in contact mechanics and turbulent flows where previous neural-physics approaches struggle. Videos and codes are available on the project page: http://www.cs.columbia.edu/cg/INSR-PDE/, Comment: ICML 2023. Project page: http://www.cs.columbia.edu/cg/INSR-PDE/

The long runtime of high-fidelity partial differential equation (PDE) solvers makes them unsuitable for time-critical applications. We propose to accelerate PDE solvers using reduced-order modeling (ROM). Whereas prior ROM approaches reduce the dimensionality of discretized vector fields, our continuous reduced-order modeling (CROM) approach builds a low-dimensional embedding of the continuous vector fields themselves, not their discretization. We represent this reduced manifold using continuously differentiable neural fields, which may train on any and all available numerical solutions of the continuous system, even when they are obtained using diverse methods or discretizations. We validate our approach on an extensive range of PDEs with training data from voxel grids, meshes, and point clouds. Compared to prior discretization-dependent ROM methods, such as linear subspace proper orthogonal decomposition (POD) and nonlinear manifold neural-network-based autoencoders, CROM features higher accuracy, lower memory consumption, dynamically adaptive resolutions, and applicability to any discretization. For equal latent space dimension, CROM exhibits 79$\times$ and 49$\times$ better accuracy, and 39$\times$ and 132$\times$ smaller memory footprint, than POD and autoencoder methods, respectively. Experiments demonstrate 109$\times$ and 89$\times$ wall-clock speedups over unreduced models on CPUs and GPUs, respectively. Videos and codes are available on the project page: https://crom-pde.github.io

This work proposes a model-reduction approach for the material point method on nonlinear manifolds. Our technique approximates the $\textit{kinematics}$ by approximating the deformation map using an implicit neural representation that restricts deformation trajectories to reside on a low-dimensional manifold. By explicitly approximating the deformation map, its spatiotemporal gradients -- in particular the deformation gradient and the velocity -- can be computed via analytical differentiation. In contrast to typical model-reduction techniques that construct a linear or nonlinear manifold to approximate the (finite number of) degrees of freedom characterizing a given spatial discretization, the use of an implicit neural representation enables the proposed method to approximate the $\textit{continuous}$ deformation map. This allows the kinematic approximation to remain agnostic to the discretization. Consequently, the technique supports dynamic discretizations -- including resolution changes -- during the course of the online reduced-order-model simulation. To generate $\textit{dynamics}$ for the generalized coordinates, we propose a family of projection techniques. At each time step, these techniques: (1) Calculate full-space kinematics at quadrature points, (2) Calculate the full-space dynamics for a subset of `sample' material points, and (3) Calculate the reduced-space dynamics by projecting the updated full-space position and velocity onto the low-dimensional manifold and tangent space, respectively. We achieve significant computational speedup via hyper-reduction that ensures all three steps execute on only a small subset of the problem's spatial domain. Large-scale numerical examples with millions of material points illustrate the method's ability to gain an order of magnitude computational-cost saving -- indeed $\textit{real-time simulations}$ -- with negligible errors.

Neural network applications have become popular in both enterprise and personal settings. Network solutions are tuned meticulously for each task, and designs that can robustly resolve queries end up in high demand. As the commercial value of accurate and performant machine learning models increases, so too does the demand to protect neural architectures as confidential investments. We explore the vulnerability of neural networks deployed as black boxes across accelerated hardware through electromagnetic side channels. We examine the magnetic flux emanating from a graphics processing unit's power cable, as acquired by a cheap $3 induction sensor, and find that this signal betrays the detailed topology and hyperparameters of a black-box neural network model. The attack acquires the magnetic signal for one query with unknown input values, but known input dimensions. The network reconstruction is possible due to the modular layer sequence in which deep neural networks are evaluated. We find that each layer component's evaluation produces an identifiable magnetic signal signature, from which layer topology, width, function type, and sequence order can be inferred using a suitably trained classifier and a joint consistency optimization based on integer programming. We study the extent to which network specifications can be recovered, and consider metrics for comparing network similarity. We demonstrate the potential accuracy of this side channel attack in recovering the details for a broad range of network architectures, including random designs. We consider applications that may exploit this novel side channel exposure, such as adversarial transfer attacks. In response, we discuss countermeasures to protect against our method and other similar snooping techniques., Comment: 14 pages, accepted to USENIX Security 2022

The biharmonic equation with Dirichlet and Neumann boundary conditions discretized using the mixed finite element method and piecewise linear (with the possible exception of boundary triangles) finite elements on triangular elements has been well-studied for domains in R2. Here we study the analogous problem on polyhedral surfaces. In particular, we provide a convergence proof of discrete solutions to the corresponding smooth solution of the biharmonic equation. We obtain convergence rates that are identical to the ones known for the planar setting. Our proof focuses on three different problems: solving the biharmonic equation on the surface, solving the biharmonic equation in a discrete space in the metric of the surface, and solving the biharmonic equation in a discrete space in the metric of the polyhedral approximation of the surface. We employ inverse discrete Laplacians to bound the error between the solutions of the two discrete problems, and generalize a flat strategy to bound the remaining error between the discrete solutions and the exact solution on the curved surface.

In this report we discuss and propose a correction to a convergence and stability issue occurring in the work of Da et al.[2015], in which they proposed a numerical model to simulate soap bubbles., Comment: 3 pages, 3 figures, technical report of Columbia University

Current quadratic smoothness energies for curved surfaces either exhibit distortions near the boundary due to zero Neumann boundary conditions, or they do not correctly account for intrinsic curvature, which leads to unnatural-looking behavior away from the boundary. This leads to an unfortunate trade-off: one can either have natural behavior in the interior, or a distortion-free result at the boundary, but not both. We introduce a generalized Hessian energy for curved surfaces, expressed in terms of the covariant one-form Dirichlet energy, the Gaussian curvature, and the exterior derivative. Energy minimizers solve the Laplace-Beltrami biharmonic equation, correctly accounting for intrinsic curvature, leading to natural-looking isolines. On the boundary, minimizers are as-linear-as-possible, which reduces the distortion of isolines at the boundary. We discretize the covariant one-form Dirichlet energy using Crouzeix-Raviart finite elements, arriving at a discrete formulation of the Hessian energy for applications on curved surfaces. We observe convergence of the discretization in our experiments., Comment: 17 pages, 18 figures

In geometry processing, smoothness energies are commonly used to model scattered data interpolation, dense data denoising, and regularization during shape optimization. The squared Laplacian energy is a popular choice of energy and has a corresponding standard implementation: squaring the discrete Laplacian matrix. For compact domains, when values along the boundary are not known in advance, this construction bakes in low-order boundary conditions. This causes the geometric shape of the boundary to strongly bias the solution. For many applications, this is undesirable. Instead, we propose using the squared Frobenious norm of the Hessian as a smoothness energy. Unlike the squared Laplacian energy, this energy's natural boundary conditions (those that best minimize the energy) correspond to meaningful high-order boundary conditions. These boundary conditions model free boundaries where the shape of the boundary should not bias the solution locally. Our analysis begins in the smooth setting and concludes with discretizations using finite-differences on 2D grids or mixed finite elements for triangle meshes. We demonstrate the core behavior of the squared Hessian as a smoothness energy for various tasks., Comment: 12 pages, 16 figures

Robotic manipulation of deformable objects is a difficult problem especially because of the complexity of the many different ways an object can deform. Searching such a high dimensional state space makes it difficult to recognize, track, and manipulate deformable objects. In this paper, we introduce a predictive, model-driven approach to address this challenge, using a pre-computed, simulated database of deformable object models. Mesh models of common deformable garments are simulated with the garments picked up in multiple different poses under gravity, and stored in a database for fast and efficient retrieval. To validate this approach, we developed a comprehensive pipeline for manipulating clothing as in a typical laundry task. First, the database is used for category and pose estimation for a garment in an arbitrary position. A fully featured 3D model of the garment is constructed in real-time and volumetric features are then used to obtain the most similar model in the database to predict the object category and pose. Second, the database can significantly benefit the manipulation of deformable objects via non-rigid registration, providing accurate correspondences between the reconstructed object model and the database models. Third, the accurate model simulation can also be used to optimize the trajectories for manipulation of deformable objects, such as the folding of garments. Extensive experimental results are shown for the tasks above using a variety of different clothing., Comment: 21 pages, 27 figures

Robotic manipulation of deformable objects remains a challenging task. One such task is to iron a piece of cloth autonomously. Given a roughly flattened cloth, the goal is to have an ironing plan that can iteratively apply a regular iron to remove all the major wrinkles by a robot. We present a novel solution to analyze the cloth surface by fusing two surface scan techniques: a curvature scan and a discontinuity scan. The curvature scan can estimate the height deviation of the cloth surface, while the discontinuity scan can effectively detect sharp surface features, such as wrinkles. We use this information to detect the regions that need to be pulled and extended before ironing, and the other regions where we want to detect wrinkles and apply ironing to remove the wrinkles. We demonstrate that our hybrid scan technique is able to capture and classify wrinkles over the surface robustly. Given detected wrinkles, we enable a robot to iron them using shape features. Experimental results show that using our wrinkle analysis algorithm, our robot is able to iron the cloth surface and effectively remove the wrinkles., Comment: 7 pages, 6 figures in IEEE International Conference on Robotics and Automation (ICRA), Stockholm, May 2016

Robotic manipulation of deformable objects remains a challenging task. One such task is folding a garment autonomously. Given start and end folding positions, what is an optimal trajectory to move the robotic arm to fold a garment? Certain trajectories will cause the garment to move, creating wrinkles, and gaps, other trajectories will fail altogether. We present a novel solution to find an optimal trajectory that avoids such problematic scenarios. The trajectory is optimized by minimizing a quadratic objective function in an off-line simulator, which includes material properties of the garment and frictional force on the table. The function measures the dissimilarity between a user folded shape and the folded garment in simulation, which is then used as an error measurement to create an optimal trajectory. We demonstrate that our two-arm robot can follow the optimized trajectories, achieving accurate and efficient manipulations of deformable objects., Comment: 8 pages, 9 figures, Proceedings of IROS 2015

We present the first triangle mesh-based technique for tracking the evolution of general three-dimensional multimaterial interfaces undergoing complex topology changes induced by deformations and collisions. Our core representation is a non-manifold triangle surface mesh with material labels assigned to each half-face to distinguish volumetric regions. We advect the vertices of the mesh in a Lagrangian manner, and employ a complete set of collision-safe mesh improvement and topological operations that track and update material labels. In particular, we develop a unified, collision-safe strategy for handling complex topological operations acting on evolving triple- and higher-valence junctions, and a flexible method to merge colliding multimaterial meshes. We demonstrate our approach with a number of challenging geometric flows, including passive advection, normal flow, and mean curvature flow., Comment: The paper has been drastically revised and submitted to a journal that requires author anonymity; while this cannot be achieved under arXiv.org policy, we have nevertheless decided to withdraw the paper to reflect our attempt to satisfy anonymity

We present a numerical model for the dynamics of thin viscous threads based on a discrete, Lagrangian formulation of the smooth equations. The model makes use of a condensed set of coordinates, called the centerline/spin representation: the kinematical constraints linking the centerline's tangent to the orientation of the material frame is used to eliminate two out of three degrees of freedom associated with rotations. Based on a description of twist inspired from discrete differential geometry and from variational principles, we build a full-fledged discrete viscous thread model, which includes in particular a discrete representation of the internal viscous stress. Consistency of the discrete model with the classical, smooth equations is established formally in the limit of a vanishing discretization length. The discrete models lends itself naturally to numerical implementation. Our numerical method is validated against reference solutions for steady coiling. The method makes it possible to simulate the unsteady behavior of thin viscous jets in a robust and efficient way, including the combined effects of inertia, stretching, bending, twisting, large rotations and surface tension.

An asynchronous, variational method for simulating elastica in complex contact and impact scenarios is developed. Asynchronous Variational Integrators (AVIs) are extended to handle contact forces by associating different time steps to forces instead of to spatial elements. By discretizing a barrier potential by an infinite sum of nested quadratic potentials, these extended AVIs are used to resolve contact while obeying momentum- and energy-conservation laws. A series of two- and three-dimensional examples illustrate the robustness and good energy behavior of the method.

Asynchronous Variational Integrators (AVIs) have demonstrated long-time good energy behavior. It was previously conjectured that this remarkable property is due to their geometric nature: they preserve a discrete multisymplectic form. Previous proofs of AVIs' multisymplecticity assume that the potentials are of an elastic type, i.e., specified by volume integration over the material domain, an assumption violated by interaction-type potentials, such as penalty forces used to model mechanical contact. We extend the proof of AVI multisymplecticity, showing that AVIs remain multisymplectic under relaxed assumptions on the type of potential. The extended theory thus accommodates the simulation of mechanical contact in elastica (such as thin shells) and multibody systems (such as granular materials) with no drift of conserved quantities (energy, momentum) over long run times, using the algorithms in [3]. We present data from a numerical experiment measuring the long time energy behavior of simulated contact, comparing the method built on multisymplectic integration of interaction potentials to recently proposed methods for thin shell contact., Comment: 10 pages, 1 figure. Accompanying report for the following publication: Asynchronous Contact Mechanics. David Harmon, Etienne Vouga, Breannan Smith, Rasmus Tamstorf, Eitan Grinspun, SIGGRAPH (ACM Transactions on Graphics), 2009

In this paper, we derive a variational integrator for certain highly oscillatory problems in mechanics. To do this, we take a new approach to the splitting of fast and slow potential forces: rather than splitting these forces at the level of the differential equations or the Hamiltonian, we split the two potentials with respect to the Lagrangian action integral. By using a different quadrature rule to approximate the contribution of each potential to the action, we arrive at a geometric integrator that is implicit in the fast force and explicit in the slow force. This can allow for significantly longer time steps to be taken (compared to standard explicit methods, such as St\"ormer/Verlet) at the cost of only a linear solve rather than a full nonlinear solve. We also analyze the stability of this method, in particular proving that it eliminates the linear resonance instabilities that can arise with explicit multiple-time-stepping methods. Next, we perform some numerical experiments, studying the behavior of this integrator for two test problems: a system of coupled linear oscillators, for which we compare against the resonance behavior of the r-RESPA method; and slow energy exchange in the Fermi--Pasta--Ulam problem, which couples fast linear oscillators with slow nonlinear oscillators. Finally, we prove that this integrator accurately preserves the slow energy exchange between the fast oscillatory components, which explains the numerical behavior observed for the Fermi--Pasta--Ulam problem., Comment: 17 pages, 3 figures, to appear in SIAM Multiscale Modeling and Simulation. v2: added new section, new figure, minor edits; v3: minor revision and clarification