1. Introduction
Clothing is ubiquitous in our daily lives. With the widespread appearance of clothing in the fashion industry, film industry, computer animation, and video games, simulating cloth has been an active research topic for more than two decades. Today, research advancement in cloth simulation has unlocked various applications such as virtual tryon (Guan et al., 2012), garment design (Umetani et al., 2011; Wang, 2018; Bartle et al., 2016; Montes et al., 2020), fold design (Li et al., 2018c), garment grading (Brouet et al., 2012), saggingfree inversion (Ly et al., 2018), and robotassisted dressing (Clegg et al., 2020; Clegg et al., 2018).
Inspired by the recent development of differentiable physics simulation and its success in rigidbody systems (de Avila BelbutePeres et al., 2018a; Geilinger et al., 2020), fluidic systems (Hu et al., 2020; Du et al., 2020), and deformablebody systems (Hahn et al., 2019; Geilinger et al., 2020)
, we argue in this paper that a large number of clothrelated applications would also benefit from a highquality differentiable cloth simulator. The critical ingredient in previous differentiable simulators is their ability to compute gradients by backpropagating any differentiable performance metrics through simulation. Such additional gradient information unlocks gradientbased continuous optimization methods, which often bring a substantial speedup compared with their gradientfree counterparts in downstream applications.
Compared with the recent active development of differentiable simulators for rigidbody and softbody dynamics, research work about differentiable cloth simulation is still relatively sparse. Indeed, cloth simulation introduces unique challenges from its codimensional dynamics and, in particular, rich contact events. While many differentiable simulators have provided solutions to derive gradients for contact models of varying complexity, their techniques typically do not expect contact to be as frequent as in cloth simulation. Deriving gradients with frequent contact and selfcollisions in cloth simulation is not fully resolved even in the stateoftheart differentiable cloth simulators (Liang et al., 2019), and our work attempts to fill this gap.
In this work, we present a differentiable cloth simulator with extra care of its contact model. We base our method on the stateoftheart cloth simulator from Ly et al. (2020) and employ its Projective Dynamics (PD) simulation method and dry frictional contact described by the SignoriniCoulomb law. Therefore, our differentiable cloth simulator inherits both the speedup from Projective Dynamics and the physical accuracy from the dry frictional contact model. Unlike previous papers relying on automatic differentiation tools to derive gradients, we present an iterative solver modified from Du et al. (2021) to accommodate the dry frictional contact model. We show that the modified iterative solver leads to a substantial speedup over a standard linear solver in gradient computation.
To have welldefined gradients, a differentiable simulator expects all quantities computed in simulation to be sufficiently smooth. However, the nonsmooth position and force changes from large numbers of contact events question this fundamental assumption. To fully understand the usefulness of a differentiable simulator in contactrich environments, our work conducts comprehensive evaluation and analysis of the behavior of our differentiable cloth simulation with a varying number of contact events. While previous papers have provided similar discussions, they primarily focus on understanding penaltybased contact models (Geilinger et al., 2020) or discussing collisions in rigidbody dynamics (Werling et al., 2021) where contact events are not as frequent as in cloth simulation. To our best knowledge, our work is the first to present such evaluation and discussion of the usefulness of gradients in a differentiable cloth simulator with dry frictional contact.
We demonstrate the efficacy of our simulator in various applications, including system identification of frictional coefficients in cloth simulation, inverse garment design for computer animation, and motion planning of robotic manipulators in robotassisted dressing. Many of these applications would either be impossible or require a much longer time to solve with existing methods. With the extra gradient information from our differentiable simulator, we unlock gradientbased optimizers to solve these problems with a much higher sample efficiency than traditional gradientfree methods.
To summarize, our work contributes the following:

We present a novel differentiable cloth simulator with dry frictional contact and an iterative solver for speeding up its gradient computation.

We evaluate the source of nondifferentiability in the dry frictional contact model and discuss the usefulness of gradients in differentiable, contactrich cloth simulation.

We show the efficacy of our simulator in various applications, including system identification, trajectory optimization for assisted dressing, closedloop control, inverse design, and realtosim transfer.
2. Related Work
Our work is closely related to cloth simulation and its applications in computer graphics. It is also relevant to the more recent differentiable simulation methods developed in the graphics and machine learning communities.
Cloth simulation
Physicsbased cloth simulation has been a popular topic in the graphics community for decades since clothing has been widely used in our daily lives. The implicit Euler integration is used to simulate cloth robustly with large time steps (Terzopoulos et al., 1987; Baraff and Witkin, 1998) while introducing excessive numerical damping. To alleviate this problem, implicit and explicit methods (IMEX) (Bridson et al., 2003; Stern and Grinspun, 2009) explicitly integrate elastic forces and implicitly integrate damping forces. Researchers have also introduced other variational integrators, e.g., BDF2 (Choi and Ko, 2002) and Symplectic integrator (Stern and Desbrun, 2006), to conserve the total energy of the system.
For cloth modeled as a massspring system, Liu et al. (2013) treat the implicit Euler integration as an energy minimization problem. The global linear system then remains constant on runtime, and each spring constraint can be solved separately in a local step. That global/local solver idea is generalized as Projective Dynamics (PD) (Bouaziz et al., 2014), which supports material models whose elastic energy has a specific quadratic form. Projective Dynamic is further extended to support more general materials (Overby et al., 2017; Liu et al., 2017). Though a local step in PD can be processed in parallel, the global step still needs to maintain a large prefactorized sparse matrix and do back substitutions in each step. Another relevant idea is PositionBased Dynamics (PBD) (Müller et al., 2007; Macklin et al., 2016), which iteratively projects each constraint in a nonlinear GaussSeidellike fashion, leading to highly parallelizable computation on GPUs.
Both PD and PBD have led to a few followup works on more advanced speedup techniques. Komaritzan and Botsch (2018) present techniques to speed up PD with contact for physicsbased character skinning. For PD and PBD in cloth simulation, Wang (2015) and a followup work (Wang and Yang, 2016) propose a Chebyshev acceleration technique that can be applied to both PD and PBD to speed up convergence. Fratarcangeli et al. (2016) also introduce a parallel randomized GaussSeidel method that reorganizes the unknowns of the sparse linear system into a few independent blocks, which can be solved in parallel with a single GaussSeidel step. Recently, the computation of integration is further accelerated by the multigrid method, e.g., geometric multigrid scheme (Wang et al., 2018), algebraic multigrid scheme (Tamstorf et al., 2015), and Galerkin multigrid scheme (Xian et al., 2019).
Lastly, our work is also relevant to previous attempts to matching cloth simulation with real fabric behaviors (Miguel et al., 2012; Wang et al., 2011; Miguel et al., 2013; Clyde et al., 2017), which is typically done by fitting constitutive material models with real material properties and can benefit from extra gradient information from a differentiable simulator (Hahn et al., 2019).
Cloth contact and friction
Contact and friction are key ingredients in modern cloth simulation. Provot (1997) proposes a method called “impact zones” to collect the nodes involved in multiple collisions into impact zones, which are treated as rigid bodies. Harmon et al. (2008) improve the failsafe of impact zones by allowing some sliding motion of the incriminated vertices. Bridson et al. (2002) introduce a hybrid method to combine the idea of applying repulsion impulses, frictions, and impact zones to handle cloth collision robustly, which is still widely used in modern cloth simulators (Narain et al., 2012; Li et al., 2020b). To model friction fully implicitly, Brown et al. (2018) treat friction as an additional dissipative term in optimization. Regarding contact handling, repulsive forces or penalty methods have been widely used (Bouaziz et al., 2014; Wang, 2015; Geilinger et al., 2020; Macklin et al., 2020) since they are easy to implement. However, these methods often need high stiffness in the penalty energy, leading to disturbing jittering artifacts that demand careful tuning. Recently, Li et al. (2020a); Li et al. (2021) define a smooth dissipative potential for friction using barrier methods. Their method can guarantee interpenetrationfree states after each time step without the need for parameter tuning.
Unlike penaltybased methods, constraintbased collision handling methods formulate contact as constraints in the physics system. Otaduy et al. (2009) first formulate cloth contacts as a sparse linear complementarity problem (LCP). Their key idea is to interleave frictional contact iterations with normal contact iterations. Instead of using a pyramidal approximation of the Coulomb friction cone (Otaduy et al., 2009), Li et al. (2018b) rely upon the exact Coulomb friction cone and adaptively refine nodes to ensure an accurate treatment of frictional contact. Although constraintbased methods typically ensure the physicsbased constraints characterizing contact and friction are satisfied after time integration, their computation cost is typically much more expensive than a penaltybased method. Recently, Ly et al. (2020) propose an efficient algorithm to incorporate frictional contact into Projective Dynamics so that nonpenetration and Coulomb constraints are satisfied simultaneously in a semiimplicit way.
Inverse dynamics
Inverse dynamics have been studied in robotics for decades to reconstruct internal forces or torques from the observations of robotic systems. However, existing methods usually focus on rigidbody systems only, which have less than a hundred degrees of freedom (DoFs) (Mistry et al., 2010; Dario Bellicoso et al., 2016; Kang et al., 2021). Inverse dynamics for highDoF systems like soft bodies, fluids, and cloth are still under exploration due to the lack of highquality numerical solutions in robotics for both simulation and differentiation. One noticeable distinction between inverse dynamics and differentiable simulation is that differentiable simulation computes additional gradients for initial states, system parameters, and design parameters. Therefore, differentiable simulation enables more applications like system identification, inverse design, and realtosim matching that traditional inverse dynamics typically do not consider.
Differentiable simulation
Differentiable simulation is a relatively recent concept explored in the graphics and machine learning community, but its original idea can be traced back to much earlier works in graphics decades ago. Perhaps one of the earliest such papers is Witkin and Kass (1988), which shows optimizing simulation to minimize an objective. Despite the recent advances in differentiable simulators in rigidbody dynamics (Popović et al., 2003; de Avila BelbutePeres et al., 2018b; Toussaint et al., 2018; Degrave et al., 2019; Qiao et al., 2020; Xu et al., 2021), softbody dynamics (Hu et al., 2019, 2020; Hahn et al., 2019; Geilinger et al., 2020; Du et al., 2021), and fluid dynamics (Treuille et al., 2003; McNamara et al., 2004; Wojtan et al., 2006; Schenck and Fox, 2018; Holl et al., 2020)
, differentiable cloth simulation still lacks a good solution. One natural idea is to use particlebased strategies that approximate a nodal system of cloth with graph neural networks
(Li et al., 2019; SanchezGonzalez et al., 2020). Although the neural networks are naturally differentiable, the physical accuracy is hard to guarantee.To accurately predict the behavior of realworld objects, recent papers (Hahn et al., 2019; Geilinger et al., 2020) present differentiable softbody systems and massspring systems with implicit Euler time integration and penaltybased contacts. Although a few recent papers (Macklin et al., 2020; Li et al., 2020a; Li et al., 2021) have introduced differentiable contact handling methods, none of them show a clothing example using the differentiability besides demonstrating the differentiability of their methods in theory. Liang et al. (2019) are the first to introduce a fully functional differentiable cloth simulation with contact, friction, and selfcollision. Instead of constructing a static LCP problem, they develop a quadratic programming (QP) problem to minimize the change between the collisionfree state and the original mesh state. Murthy et al. (2021) also present a differentiable cloth simulator with penaltybased frictional contact in their fully differentiable simulation and rendering pipeline. In this paper, we build upon the differentiable PD framework (Du et al., 2021) to simulate cloth dynamics with dry frictional contact (Ly et al., 2020) and augment it with gradient computation.
A significant challenge in differentiable simulation is the presence of nonsmooth contact events, where nonsmoothness can arise from discontinuous contact shapes (Popović et al., 2000), impulsive forces (Hu et al., 2019), or branches in contact laws (Ly et al., 2020), to name a few. Many existing differentiable simulators assume such nonsmoothness does not affect the usability of gradients. Indeed, they present empirical results that observe the benefits of using gradients in downstream applications (de Avila BelbutePeres et al., 2018a; Du et al., 2021; Liang et al., 2019; Geilinger et al., 2020)
. While our work also observes the advantages of gradients over gradientfree methods in several applications, we take one step further by analyzing the source of nonsmoothness and discontinuities in contactrich cloth simulation. We hope that our experience can serve as useful heuristics for applying differentiable simulation methods in the future.
3. Background
In this section, we briefly review the Projective Dynamics cloth simulation method with dry frictional contact described in Ly et al. (2020), on which our differentiable cloth simulator is based. Our core contribution lies in the development of gradient computation in this PD framework with dry frictional contact, which we detail in the next section.
3.1. Projective Dynamics Method
Implicit time integration
We model cloth as a nodal system with 3D nodes. Let and
be two vector functions from
to indicating the positions and velocities of all nodes at time . We consider the standard implicit timestepping scheme discretizing and as follows:(1)  
(2) 
where is the time step size, a positive diagonal mass matrix of size , and and the internal and external forces at each node stacked as two dimensional vectors, respectively. Note that we assume for now and do not contain contact forces, which will be discussed separately in Sec. 3.2. We use the cloth material model described in Bouaziz et al. (2014) to define . The implicit timestepping scheme connects the state of the nodal system at time to at new time .
Optimization view
As discussed by Martin et al. (2011), the implicit timestepping scheme can be rephrased as finding a saddle point of the following energy minimization problem:
(3) 
where is a constant vector that can be precomputed at the beginning of the current time step. The potential energy defines the internal force by its spatial gradients: . It is easy to verify that setting the gradient of the objective function to zero leads to a system of equations identical to Eqns. (1) and (2).
Local and global solvers in PD
The key assumptions in PD is that the internal energy can be written as a sum of quadratic forms (Bouaziz et al., 2014; Liu et al., 2017):
(4)  
(5) 
Here, each energy projects
, a linear transformation of
, to its closest point in the set and scales its squared distance by a prespecified stiffness . Both and are predefined and independent of . Here, is the constraint manifold, and is the auxiliary projection variable as defined in Bouaziz et al. (2014).With such an assumption on , PD proposes a localglobal solver to minimize a surrogate objective function:
(6) 
where stacks up all from each . In the local step, PD fixes and projects each to its corresponding by solving Eqn. (4). Such a local step can be done in parallel for each . In the global step, PD fixes and minimizes as a function of , which turns out to have a closedform solution:
(7) 
By alternating between the local and global steps, PD monotonically decreases the surrogate energy until convergence, which can be shown to agree with the saddle point of (Liu et al., 2017). The source of efficiency in PD comes from the observation that is a constant matrix that can be prefactorized at the beginning of the simulation, leading to an efficient global step requiring backsubstitution only.
3.2. Dry Frictional Contact in Projective Dynamics
SignoriniCoulomb law
Ly et al. (2020) augment the standard PD framework described above with nonpenetration collision and Coulomb friction by assuming contact applies to nodes only. At each time step, assuming a collision detection algorithm has identified a set which describes the indices of contact nodes. For each node , the SignoriniCoulomb law (Brogliato, 2016) requires its local force and velocity to satisfy one of the following three conditions:
(8a)  
(8b)  
(8c) 
Here, is the frictional coefficient, and and are the nodal velocity and contact force represented in the local contact frame spanned by the tangential plane and contact normal on the contact surface. The notations and represent the tangential and normal components of a 3D vector defined in the local frame. We further define as the set of valid pairs described by the conditions above, allowing us to rewrite Eqn. 8 compactly: .
Implicit time integration with dry frictional contact
The original implicit time integration now needs to be augmented with additional constraints describing contact conditions (Ly et al., 2020):
(9) 
Here, we have added the subscript in the definitions of the contact set and the contact condition described above to specify the time step they are defined from. The notation and select the 3D force and velocity corresponding to the th node from and , respectively. The Jacobian matrix of size selects global vectors defined on the contact nodes and computes their coordinates in the local contact frames. Note that our definition of implies that it is computed in an explicit manner, i.e., it relies on the last state entering the th time step. For brevity, this background section will describe simple contact events between a node and a static plane, in which case is a constant matrix that does not need to be computed from . Extensions to more sophisticated contact events, e.g., selfcollisions between multiple contact nodes or contact events between a node and a nonplanar obstacle, can be found in Ly et al. (2020) and in our implementation and experiments.
Projective Dynamics with dry frictional contact
The core idea in Ly et al. (2020) is an additional local step in PD that solves each contact node independently. Noting that the contact conditions are primarily defined on nodal velocities instead of nodal positions, Ly et al. (2020) first rewrite the global step in Eqn. (7) as an equation of velocities:
(10) 
The velocitybased PD globallocal steps then alternate between this new global solver and the original local step, which is unaffected by this change of variables. By comparing Eqn. (9) and Eqn. (10), we can see that incorporating the contact force adds an additional impulse to the righthand side of the global step:
(11) 
The goal of the global solver is now to find that satisfy contact conditions at each . Noting that in Eqn. (11), is the only operator that couples unknown contact forces from different , Ly et al. (2020) propose the following iterative solver based on the decomposition of to solve Eqn. (11):
(12) 
where the superscripts and indicate two consecutive iterations in the iterative solver. It is easy to see that the proposed iterative solver fully decouples the momentum equation on each node, allowing Ly et al. (2020) to enforce the SignoriniColumb law by adjusting on each contact node independently in a straightforward manner. When the iterative solver converges, Ly et al. (2020) prove that the result is indeed a solution to Eqn. (9).
4. Differentiable Cloth Simulation
We now describe in detail how we extend the forward simulator in Sec. 3 to build a differentiable cloth simulator. We start by differentiating implicit time integration in PD without contact, followed by explaining how contact gradients can be added to this framework. Compared with other differentiable simulators, our simulator is unique because of its treatment of the contactrich nature of cloth simulation: many existing differentiable simulators focus on sparse contact events in rigidbody or deformablebody dynamics (Geilinger et al., 2020; Hu et al., 2019; Du et al., 2021), and the stateoftheart differentiable cloth simulation method (Liang et al., 2019) handles rich contact events but without physicsbased contact forces or friction. To our best knowledge, our work is the first to present a differentiable cloth simulator that can handle rich contact events with Coulomb’s law of friction.
4.1. Gradients without Contact
Differentiating implicit time integration
The core step in building a differentiable simulator based on implicit timestepping scheme is to backpropgate gradients through the implicit integration described in Eqns. (1) and (2), or equivalently, to derive the Jacobian of the output with respect to the input . Mature techniques such as sensitivity analysis, adjoint method, and implicit function theorem have proven to be successful in computing such gradients (Hahn et al., 2019; Geilinger et al., 2020; Du et al., 2021). To sketch the idea, we plug from Eqn. (2) into Eqn. (1):
(13)  
which is essentially the firstorder optimality condition in Eqn. (3). Differentiating on both sides gives
(14)  
or equivalently,
(15)  
In backpropagation, such a Jacobian matrix is coupled with the gradients of a loss function
which are passed to the previous state (assuming the gradients below are both column vectors):(16)  
Here, we obtain the adjoint vector by solving the linear system of equations with the matrix on the lefthand side, which avoids an explicit inversion of the large and sparse matrix. Backpropagating gradients of with respect to and can be derived similarly. In fact, it solves the same linear system but with a different vector on the righthand side.
Differentiating with Projective Dynamics
With assumptions on the energy form in PD, Du et al. (2021) show that we can speed up the computation in by exploiting the special form of . The first and secondorder derivatives in Eqn. (5) are given by the following equations (Du et al., 2021; Liu et al., 2017):
(17)  
(18) 
where is the projection of onto obtained in the local step. Interested readers can refer to the appendix of Liu et al. (2017) for their derivation details. We plug them to and rewrite the linear system as follows:
(19) 
Noting that a prefactorization of exists in Projective Dynamics, Du et al. (2021) propose the following iterations to solve :
(20) 
Similar to the localglobal steps in Projective Dynamics, a local step can evaluate across each in parallel, and a global step, which reuses , can solve with backsubstitution only. Du et al. (2021) show that such a localglobal solver is empirically faster than directly solving Eqn. (16), which resembles the source of efficiency in the original PD method.
4.2. Gradients with Contact
While Du et al. (2021) have discussed extensively how to fuse backpropagation into the Projective Dynamics framework, its support of contact is limited to nonpenetration conditions without a proper treatment on friction. On the other hand, other prior papers have provided solutions to deriving gradients with contact governed by SignoriniCoulomb law, but they focus on either smallscale problems (e.g., rigid bodies in de Avila BelbutePeres et al. (2018a)) or sparse contact events on a static ground only (Geilinger et al., 2020). Here, we present our differentiable cloth simulator that both inherits the speedup from PD and handles gradients with complicated contact events like selfcollisions, which are overlooked in many existing differentiable simulators but common in cloth simulation.
Differentiating implicit time integration with contact
Consider the th time step in simulation with contact which takes as input and computes that satisfy Eqn. (9). In particular, for each contact node , the forward simulation identifies which one of the contact conditions in Eqn. (8) applies to . As an example, assuming a contact node satisfies Eqn. (8c), its constraints can be summarized as follows:
(21a)  
(21b)  
(21c)  
(21d) 
where and extract the and components of a twodimensional vector, respectively. The other two cases of Eqn. (8) similarly enforce three equality constraints ( for taking off and for sticking) and a number of inequality constraints on . If the inequality constraint is inactive, slightly perturbing the inputs to the simulator will keep
inside its interior. Therefore, we can remove the inactive inequality when deducing the gradients during backpropagation. When the inequality constraint is active, the gradients of the simulation are not well defined because it represents corner cases that can be categorized into more than one contact types. These corner cases introduce nonsmoothness to the simulator, but it is worth mentioning that they do not create discontinuities, just like the standard rectifier (ReLU) activation function is still continuous despite the nonsmoothness at its turning point.
After removing the inequality constraint, we further define a nonlinear vector function for the lefthand side of the three equality constraints and compactly represent the constraints as . It is convenient that for all three cases in Eqn. (8) have three dimensional outputs. We can stack from all contact nodes into a nonlinear function , allowing us to restate Eqn. (9) as follows:
(22) 
Note that we choose instead of as our variable in order to be consistent with the forward PD simulation with contact described by Ly et al. (2020). Differentiating both sides of the first equation with respect to leads to the following result:
(23)  
(24)  
(25) 
Comparing it with Eqn. (15), we see the matrix to be inverted now has an additional component dependent on , which we obtain from differentiating the constraint in Eqn. (22):
(26)  
(27) 
We stress that computing and is trivial because both partial derivatives are blockdiagonal matrices. Therefore, can be parallelized among all contact nodes. Backpropagation through to can be implemented with the same adjoint method before, which we give in the equation below for completeness with the notation overloaded:
(28) 
Speedup with Projective Dynamics
With additional information about from PD, we can rewrite the adjoint vector in Eqn. (28) by comparing it with Eqn. (19):
(29) 
which naturally leads to the following iterative iteration:
(30) 
Comparing it with Eqn. (20), we see the role of is replaced with . Since we have shown that computing can be massively parallelized among contact nodes, it is suitable for contactrich cloth simulation and preserves the efficiency of the localglobal solver.
Convergence
In theory, such an iterative solver is guaranteed to converge from any initial guesses of when the spectral radius . Empirically, we notice in our cloth simulation that divergence is uncommon, especially when highprecision backpropagation is not required (Sec. 6). When the iterative solver fails to converge, we switch back to a direct sparse matrix solver to solve Eqn. (29).
Extensions
We deliberately skipped the full definition of in all equations above for a clearer presentation of the main idea behind our differentiable cloth simulation. We now elaborate on the role of and discuss its implications on more complex collisions. When the contact surface is static but nonplanar with spatially varying surface normals,
is dependent on the positions of the nodes where the contact events occur. In this case, we estimate the contact normal based on the positions where the contact events occur, which are given by any collision detection algorithms. Replacing
with requires very minimal extra work to the gradients in Sec. 4.2, as it contributes to the gradients with respect to in a straightforward manner.Another type of collisions we consider is selfcollisions between nodes. In this case, contact occurs between pairs of nodes where the contact normals are defined by the relative position between two nodes in the pair. Therefore, we can still define as a function of with each row block now consisting of two blocks of nonzero elements corresponding to the two contact nodes in the pair. The gradient derivation in Sec. 4.2 remains unchanged in this case except that the dependencies between and need to be added to Eqn. (23). As we inherit from Ly et al. (2020) the contact model on nodes only, we also inherit its limitation of not handling other types of selfcollisions like edgeedge collisions or vertexface collisions. We leave the derivation of gradients with such cases as future work.
5. Evaluations
This section evaluates and discusses numerical properties of the proposed differentiable cloth simulation method in Sec. 4. We first evaluate its gradients by analyzing their source of nonsmoothness and studying their usefulness in highdimensional optimization problems. Next, we compare the dry frictional model with the contact model in the stateoftheart differentiable cloth simulation method (Liang et al., 2019) and discuss their differences. Finally, we analyze the numerical properties of our iterative solver in backpropagation and compare its performance with a direct linear solver.
5.1. Continuity and Smoothness
A fundamental assumption in any differentiable simulators is that the underlying physics system is smooth so that gradients can be welldefined. Eqn. (9) contains three possible sources of discontinuities and nonsmoothness in the order of their occurrences: determining the contact set , computing that represents the local contact frames, and choosing between three branches from in Eqn. (8). Below we discuss the effects from each of them in detail.
Continuity of branches in contact conditions
First, we hypothesize that switching between the three branches in Eqn. (8) does not create discontinuity in Eqn. (9): consider the th time step and assume that is fixed and that is a continuous function of , if we treat Eqn. (9) as a function that takes as input and returns , we hypothesize this function is continuous. In other words, perturbing a little will not cause jumps in the resultant even though the corresponding may need to switch between branches during the perturbation. The intuition is that the three branches in Eqn. (8) together define a connect set in which from one branch can smoothly transition to another.
While we hypothesize that these branches do not create discontinuities, they do introduce nonsmoothness due to corner cases that can be arbitrarily classified into either branch, e.g., when a node is static but about to slip. Gradients at these corner cases are not welldefined, but the subset these corner cases occupy in
has measure zero. Therefore, we still expect gradientbased optimization to be functional, just like we have observed the success of gradientbased optimization in modern deep learning with nonsmooth but continuous operators, e.g., ReLU, max pooling, and so on.
We empirically validate this hypothesize with the following experiment. We simulate a piece of cloth on a rigid, static, and frictional sphere for 200 timesteps. The sphere has one frictional coefficient , and we repeat the experiments by varying from 0 to 0.35 (Fig. 2). When is large, all nodes are fixed on the sphere due to their large static friction. When is close to , the cloth slides on the sphere under gravity and takes off from the sphere near the end of the simulation. As changes from to , each node on the cloth undergoes the transition from having static friction (sticking) to having dynamic friction (slipping), and eventually having no friction (taking off). However, since each node has a different contact normal, the turning point for each of them to switch between these branches is different. Overall, when we gradually change , the ratio among the numbers of nodes with static friction, dynamic friction, and no friction at the end of the simulation also changes gradually, allowing us to observe how switching between these branches affects the continuity of the physical quantities in simulation.
We summarize the quantitative results from this simulation experiment in Fig. 2 (green curves). Specifically, we plot the velocity of three nodes , , and as a function of at an intermediate time step (50) and near the end of simulation (200). We select three nodes from the corners, edge centers, and the center of the cloth, respectively. All velocities converge to when becomes large, which is as expected because a large leads to static friction that freezes these nodes. We also notice a turning point in the velocity curve for each node, indicating a switch between static and dynamic friction. The turning points are located at slightly different values for each node because their normals on the sphere surface are different. We observe from the figure that the branches in Eqn. (8) do not cause discontinuities, although they may introduce nonsmoothness. This observation confirms empirically our hypothesis on discontinuities and nonsmoothness before.
Continuity of local frames at contact nodes
A second source of possible discontinuity and nonsmoothness comes from computing , which can be decomposed into two steps: determining a contact point on the contact surfaces and calculating the local normal and tangent vectors. Both steps depend on the geometric representation of the contact surfaces. As pointed out by previous papers in rigid body dynamics (Popović et al., 2000; Werling et al., 2021), contact surface discretization causes jumps in surface normals, and therefore they create discontinuities in velocity and position calculation. To confirm this observation in cloth simulation, we repeat the previous experiment by replacing the analytically described sphere with a triangulated one (pink triangulated sphere surface in Fig. 2 right). After such replacement, we clearly observe jumps in the intermediate and final velocity from the selected contact nodes (pink curves in Fig. 2). Note that the jumps in velocities from the final velocities are less evident than from the intermediate velocities because of their vertical axes have different ranges. This observation agrees with similar experiments from previous papers about rigid body dynamics and suggests that one should favor analytical surfaces in differentiable cloth simulation whenever possible.
We end our discussion on the discontinuities from with two more remarks. First, we notice the contact node velocity curves in Fig. 2 are partitioned into a small number of continuous segments. While this is consistent with the result reported by previous work (Popović et al., 2000) discussing contact events on rigid body dynamics, we still find it surprising because cloth simulation involves a larger number of degrees of freedom and our experiment involves many more jumps in contact normals from many contact events. Therefore, we expected the velocities would contain much more discontinuous points. Second, if the scene contains multiple objects where the cloth can be in contact with, it is not uncommon to see jumps in the locations of contact events, e.g., from one object to another, even with a continuous representation of each object. Such jumps naturally lead to large discontinuities in simulation. While we do not provide solutions to it in this work, a closely related issue has been extensively studied in differentiable rendering (Li et al., 2018a) from which we may borrow inspirations in the future.
Continuity of contact sets
The last and most common source of discontinuity comes from deciding the contact set at each time step, i.e., whether a node should be added to the contact set or not. Obviously, this selection process is not continuous. It is worth noting that having a constant contact set does not cause too much trouble for gradient computation regardless of its size (Fig. 2). Instead, it is the change in the contact set from time to time that brings in discontinuities because whenever a new node is put in the contact set, it adds an impulsive force to the righthand side of Eqn. (9).
To better understand the effects of changes in contact sets, we hang a piece of cloth above a static sphere in simulation and let the lower half of the cloth fall and slide on the sphere due to gravity (Fig. 3 top). We equip this experiment with a system identification task of the frictional coefficient and the stiffness parameter of the cloth: given a motion sequence of the cloth generated from a pair of unknown and , we define a loss function that sums over all time steps the squared distance between each node position in simulation and its corresponding location in the given motion sequence. We repeat the task with four settings of the sphere at different horizontal offsets, leading to a varying frequency of contact events among them.
We plot the landscape of the loss function in Fig. 3 (middle) and compare its smoothness among the four settings with different frequencies of contact events. At first glance, it seems that all four landscapes are equally smooth. However, magnifying a small region of each landscape shows that there is a profound distinction between their continuity and smoothness (Fig. 3 bottom): as establishing and breaking contact becomes more frequent, the local landscape tends to be bumpier.
Summary
We have discussed the three sources of potential discontinuities and nonsmoothness in our differentiable cloth simulator ordered by their damage to the gradients: the branches in the contact conditions only introduce nonsmoothness and maintain continuity; contact surface discretization creates discontinuities due to the jumps of normals across adjacent triangles, but we still expect continuity almost everywhere; adding or deleting nodes in the contact set creates frequent and the most severe discontinuities due to the introduction or removal of impulsive forces.
5.2. Usefulness of Gradients
One advantage of gradientbased approaches over their gradientfree counterparts is their faster convergence rate: typically, gradientfree methods explore the local landscape of an objective function by evaluating it with massive samples in the neighborhood. However, such a sampling approach quickly becomes less efficient as the dimension of the decision variables grows higher. Due to this reason, we hypothesize that using gradients from our differentiable simulator will become most beneficial when we have a large number of decision variables. We verify this hypothesis using a control optimization problem: we simulate a piece of cloth with timeinvariant external forces applied to each node. The goal is to design the force at each node to pull the center of the cloth to a target position in the end (Fig. 4 top). We simulate the cloth with four settings on its DoFs: nodes, nodes, nodes, and nodes. Therefore, the number of variables in the optimization problem is , , , and , respectively.
Fig. 4 reports in this problem the convergence rates of LBFGSB (Liu and Nocedal, 1989) and CMAES (Hansen, 2006), two representative gradientbased and gradientfree approaches, with varying DoFs. We start with the largest setting of the cloth (300 DoFs) and run both LBFGSB and CMAES with five randomly chosen initial guesses on the force values in parallel. We then plot the loss vs. time step curves for both methods in Fig. 4 bottom left. Here, the loss is the minimum loss from each method’s five parallel runs as a function of the time steps they have consumed during the optimization process. It is clear to see the obvious speedup of LBFGSB over CMAES in this 300DoF optimization problem, just as we expected.
Additionally, we repeat the experiment with three other settings of the cloth and plot the speedup vs. DoF curves in Fig. 4 bottom right, where the speedup is the ratio between the number of time steps used by CMAES and LBFGSB when their losses reach 0.01. Again, we observe significant speedup from LBFGSB across the board, with the largest speedup from the cloth with the highest DoFs. This observation confirms the benefits of using gradients from our differentiable cloth simulator in inverse design problems.
5.3. Benefits of Dry Frictional Contact
We compare the contact model in our differentiable simulator with that in the stateoftheart differentiable cloth simulator from Liang et al. (2019). Both models ensure penetrationfree simulation and detect selfcollisions (nodenode collisions in ours and vertexface and edgeedge collisions in Liang et al. (2019)). However, the contact model in Liang et al. (2019) does not take into account complementarity conditions on either contact forces or frictional forces, which may lead to undesirable artifacts. To visualize this issue, we consider a test scenario in which a napkin falls freely onto the inner surface of a bowl (Fig. 5). The differentiable simulators from both Liang et al. (2019) and our work manage to simulate the napkin without numerical explosion. However, the dry frictional contact model in our differentiable simulator shows more physically realistic motion (Fig. 5 middle), whereas the contact model from Liang et al. (2019) leads to more drastic changes in the size of the napkin with a popping artifact after its contact with the bowl (Fig. 5 top and bottom). These artifacts are explainable because the collision handling algorithm in Liang et al. (2019) modifies node positions after penetration without verifying whether such an update requires sticky contact forces. When the napkin becomes in contact with the concave inner surface of the bowl, such a direct modification injects extra elastic energies. This experiment confirms supporting a more advanced contact and friction model in differentiable cloth simulation is both viable and beneficial.
5.4. Evaluation of Iterative Solver in Backpropagation
Another key component in our differentiable cloth simulation is the iterative solver in Eqn. (29) that utilizes the prefactorized in PD. In a standard differentiable simulator, solving Eqn. (29) would be done with a sparse matrix solver treating as a whole. To understand the performance of the proposed iterative solver, we design two benchmark tests: a “Wind” test where a hanging napkin moves under synthetic wind and a “Slope” test where a ribbon slides on a slanted plane (Fig. 6). These two tests represent two extreme cases in contact handling: the napkin in “Wind” barely has any contact or selfcollisions, whereas every single node of the ribbon in “Slope” is in contact with the plane. We then compare the time cost of our iterative solver with a sparse LU solver in both tests. We implement both solvers using Eigen (Guennebaud et al., 2010) and choose LU factorization because is usually not a symmetric or positive definite matrix, preventing us from using more specialized sparse matrix solvers like Cholesky factorization or Conjugate Gradient methods. For the iterative solver, we report two results from lowprecision (1e4) and highprecision (1e6) convergence thresholds that control the termination condition in the iterations. We find these thresholds by varying it from 1e1 to 1e9 until the gradients computed from the iterative solver start to agree with the direct solver, which we treat as the groundtruth gradients. We repeat the experiments with three mesh resolutions to test the scalability of our method.
We summarize the statistics from both tests in Table 1. Overall, the iterative solver in backpropagation is faster than the direct solver, and the speedup becomes more substantial as we increase the mesh resolution. The speedup is also more evident with lowprecision threshold, which is consistent with the results reported in previous PD papers (Bouaziz et al., 2014; Liu et al., 2017; Du et al., 2021). A further decomposition of time confirms that solving the linear system takes up the majority of backpropagation time, so any improvement in the choice of solver has a dominating positive effect. Although the lowprecision results may not match the groundtruth gradients as closely as their highprecision counterparts, their backpropagation is significantly faster and can still benefit gradientbased optimization. This is because even an imperfect gradient can still guide gradient descent algorithms to converge as long as it is along the descending direction. This is confirmed empirically by our experiments in Sec. 6, in most of which we use lowprecision convergence threshold and still observe successful gradientbased optimization results. It is also worth mentioning that contactrelated operators, whose time cost is included in the “Other” and “” columns in Table 1, add very little extra overhead to the backpropagation. Finally, we observe uncommon but nonnegligible failure cases of the iterative solvers in one experiment ( with 1e6 as the threshold in Table 1) “Wind”. For the failed timesteps, we switch to the sparse LU solver, adding an extra time in backpropagation.
Test  Res.  Solver  Backprop Time (s)  (%)  (%)  Iter. (%)  Direct (%)  Other (%)  Fail (%)  Speedup 

Wind  12x12  Direct  0.606  13.9  0    84.4  1.7  0   
Ours (1e4)  0.212  39.9  55.5    4.6  2.9x  
Ours (1e6)  0.424  19.2  78.5    2.3  1.4x  
24x24  Direct  3.360  9.9  0    89.1  1.0  0    
Ours (1e4)  0.591  58.2  36.3    5.5  5.7x  
Ours (1e6)  2.858  11.7  87.2    1.1  1.2x  
48x48  Direct  24.001  6.8  0    92.5  0.7  0    
Ours (1e4)  2.262  63.2  29.7    7.1  0  10.6x  
Ours (1e6)  29.255  5.4  82.6  11.5  0.6  15  0.8x  
Slope  12x12  Direct  0.978  13.0  3.4    82.4  1.2  0   
Ours (1e4)  0.312  70.1  9.4  16.1    4.5  3.1x  
Ours (1e6)  0.643  18.9  2.7  76.8    1.6  1.5x  
24x24  Direct  5.331  10.4  1.5    87.4  0.7  0    
Ours (1e4)  0.781  70.1  9.4  16.1    4.5  6.8x  
Ours (1e6)  1.507  34.7  4.9  58.1    2.4  3.7x  
48x48  Direct  37.296  6.5  0.7    92.5  0.4  0    
Ours (1e4)  3.095  68.1  8.2  19.8    3.8  12.0x  
Ours (1e6)  6.825  31.2  3.5  63.6    1.8  5.5x 
6. Applications
In this section, we demonstrate a variety of clothrelated applications that can benefit from our proposed differentiable cloth simulation with dry frictional contact. We run all optimizations on a workstation of 80 CPU cores and 80G memory. Depending on the problem complexity, the wallclock time of running these applications varies from less than 30 minutes to 2 hours for all methods. We write the backbone of our simulator in C++ and use Eigen for matrix and vector operations. Each application defines an optimization problem which we solve with LBFGSB, a classic gradientbased optimizer that can leverage the differentiability of our cloth simulator while limiting parameters to physicallyplausible ranges. We use the implementation of LBFGS++ (Yixuan, 2021) for LBFGSB, which implements the MoréThuente line search (Moré and Thuente, 1994). We repeat all experiments with multiple random seeds and use the same random seed set (therefore same initial parameter set) for all methods. We report the optimization results and comparisons with other gradientfree algorithms in Table 2. In Table 3, We report the time steps used by each optimization method to reach minimum final loss as well as the convergence ratio of our iterative solver during backpropagation. In Table 2 and all loss vs. time steps plot below, we report the minimum loss achieved across all random seeds. Below we describe each application and highlight major results. More information of each example is detailed in Appendix A.
Optimized Loss Percentage  

LBFGSB (Ours)  CMAES  (1+1)ES  
Name  Param #  Seed #  Min Init. Loss  Final (%)  Final (%)  Equal Time (%)  Final (%)  Equal Time (%) 
Tshirt  6  5  6.62  0.53  0.60  15.78  1.74  2.43 
Sphere  1  10  0.46  0.43  0.00  0.43  0.65  2.37 
Hat  18  5  21.83  0.42  0.55  37.51  0.51  5.87 
Sock  36  5  17.08  11.59  17.39  64.70  18.43  52.80 
Dress  2  16  0.90  76.91  79.02  80.02  77.25  82.57 
Flag  8  10  2.88  4.77  38.89  10.57  40.14  25.57 
Name  Convergence Time Steps  Iter. Conv. %  

LBFGSB  CMAES  (1+1)ES  
Tshirt  4500  53500  47750  99.75 
Sphere  2450  18900  11200  100.00 
Hat  11200  159200  106000  90.58 
Sock  17600  176800  91600  99.41 
Dress  2750  3125  10000  84.80 
Flag  6800  17300  24800  96.00 
6.1. System Identification
Tshirt
In the “Tshirt” example, we are given a sequence of motions of a hanging Tshirt under synthetic wind generated from a parameterized sinusoidal function. The goal is to estimate a material parameter in the cloth (1 DoF) and identify the wind model parameters (5 DoFs controlling the amplitude, phase, and frequency of the sinusoidal function in three dimensions) from the motion data. We define the loss function as the L2distance between the nodal positions of the Tshirt from the simulation and the given motion sequence.
Unlike the “Wind” example in Sec. 5, contacts and frictions are much more frequent in this example due to selfcollision between the front and back layer of the Tshirt under the wind. We show the simulation of the Tshirt with parameters before and after optimization in Fig. 7 and compare the three optimizers in Table 2 and Table 3: LBFGSB, CMAES, and (1+1)ES. Both CMAES and (1+1)ES are standard gradientfree evolutionary strategies (ES). We can conclude from Fig. 7 and Table 2 and 3 that all three methods manage to optimize system parameters leading to motion sequences visually identical to the given input, but LBFGSB converges much faster due to the extra knowledge of gradients.
Sphere
To highlight the effect of dry frictional contact in our simulator, we create the “Sphere” example with the goal of matching the motion sequence of a cloth interacting with a sphere by estimating the frictional coefficient of the sphere. In this example, we let a piece of square cloth fall freely on a sphere, whose motion after being in contact with the sphere is largely controlled by the frictional coefficient. This example involves both selfcollisions between nodes on the cloth and external contacts with the sphere. Similar to the example above, we run both LBFGSB and two ES baselines and report their statistics in Table 2 and 3. All methods can optimize to a frictional coefficient that generates a motion sequence visually identical to the given input. However, this optimization problem is special in that it only has one parameter and the landscape of the loss function is bumpy due to the large variation in the collision set as a function of the frictional coefficient, as suggested by the Continuity of contact sets experiment in Sec. 5.1. In this case, evolutionary strategies can achieve a lower final loss because the samples can freely explore the full range of frictional contact, while LBFGSB can get trapped in a local optimum.
6.2. RobotAssisted Dressing
Another line of research that can benefit from a differentiable cloth simulator is robotassisted dressing. The mainstream solution to these tasks is typically gradientfree methods like evolutionary strategies, reinforcement learning, or inverse dynamics before a differentiable cloth simulator becomes available. We present two examples to demonstrate the usage of gradients in robotassisted dressing: “Hat” (Fig.
9) and “Sock” (Fig. 10). In both examples, the goal is to find trajectories for a kinematic robotic manipulator to put on the hat or the sock. The end effectors of the manipulator pick a few prespecified vertices on the cloth meshes and pull them along the kinematic trajectories, which are parametrized as Bsplines. By optimizing the parameters of the Bsplines (18 DoFs in “Hat” and 36 DoFs in “Sock”), we can direct the manipulator to move the hat until it reaches the target location on top of the sphere. We define the loss function in “Hat” as the L2distance between the hat’s final position at the last time step and a predefined target position, which we generate by translating the hat’s rest shape onto the top of the sphere. The loss function for “Sock” is defined as the L2distance between the desired and simulated locations of a few key points manually chosen on the sock and evaluated at the middle and the end of the simulation.With the gradient information at hand, we run the LBFGSB optimizer to tune the parameters of the trajectories and compare its performance with CMAES and (1+1)ES (Table 2 and 3). We notice that within the same time budget, LBFGSB converges substantially faster than the gradientfree baselines to a better solution (Fig. 9 and Fig. 10). We can safely conclude that the fast convergence of LBFGSB unlocked by our differentiable cloth simulator is a clear advantage over gradientfree methods.
6.3. Inverse Design
The next application we present is “Dress”, an inverse design example that aims to optimize cloth material parameters in a dress so that its dynamic motion can satisfy certain design intents. Specifically, we optimize the material parameters of a twirl dress so that after the dress spins, the apex angle of the conelike dress agrees with the target value (100 degrees in our experiment). We define the loss function as the difference between the hemline height corresponding to the target apex angle and the estimated apex angle from points on the hemline of the dress at the last frame of the simulation. We report the optimization results from LBFGSB and two ES baselines in Table 2 and visualize the simulation results before and after optimization in Fig. 1. Similar to the previous tasks, We notice that LBFGSB achieves better optimized results using fewer time steps.
6.4. A RealtoSim Example
In this section, we present a realtosim “Flag” example (Fig. 11). In this example, we use the realworld motion sequence captured on a flag flapping in the wind from previous work (White et al., 2007) and aim to reconstruct a digital twin of the scene in simulation. This includes not only estimating the material parameters of the flag but also modeling the unknown wind condition at the capture time, which is particularly challenging due to its intricate stochastic model with unknown degrees of freedom. We model the wind force at each time step as a 3D force applied near the center of the scene and spatially decaying proportional to the inverse distance to the center. To model the transient nature of the wind force, we modulate magnitude of the 3D vector of a sinusoidal wave as a function of time with parameterized frequency and phase offset. Together, the material and wind model define an eightdimensional parameter space to be optimized. We define the loss function as the L2distance between the positions of all nodes at each time step in simulation and the groundtruth motion sequence.
We solve this task using LBFGSB and the two ES methods and report their performance in Table 2, 3, and Fig. 11. All three methods substantially reduce the loss after optimization, but LBFGSB achieves a lower final loss. We plot the trajectories of 6 key nodes before and after optimization (orange) along with the groundtruth reference trajectories (yellow) from the motioncaptured data in Fig. 11. By comparing the left and right images in Fig. 11, we can see that the LBFGSB optimizer reduces the discrepancy substantially between the simulated and actual trajectories after optimization. The realtosim matching is still imperfect as indicated by the nonzero final loss, which we suspect could be due to the simplistic nature of our synthetic wind model. A more sophisticated and expressive wind model, e.g., a neural network, may serve a better role in modeling and matching the realworld physics, which we leave as future work.
6.5. Hat Controller
We end this section with an advanced “Hat” task. Unlike the previous openloop trajectory optimization with a fixed starting position, the goal of this new task is to train a generalizable closedloop controller that can put on the hat from a random starting position sampled from a fixedradius hemisphere around the head. Specifically, we train a closedloop control policy , which takes as input the current state of the task and outputs an action vector at each time step . The state vector includes the hat node positions, the orientation of the hat, and the distance between the two end effectors, and the action vector represents the position of the two end effectors at the next time step. We represent the control policy as a neural network parametrized by
consisting of two fullyconnected hidden layers with 64 neurons and
tanh activation functions (117126 parameters in total). To train the controller with gradientbased optimization, we integrate the neural network policy with our differentiable simulator as described in Fig. 13.In each epoch during training, we randomly sample 20 starting positions of the hat and compute a loss averaged from all simulation sequences. The loss function of each sequence is defined by
, where measures the stretching of the hat using the distance change between the two endeffectors, measures the L2distance between the lasttimestep pose of the hat and the target pose, and is the orientation difference between the lasttimestep pose and the target pose of the hat. The gradient of the loss is then computed by our differentiable simulation framework and used by a gradientbased solver, e.g., Adam (Kingma and Ba, 2017), to update the policy parameters. For testing, we evaluate the controller from 20 fixed starting position configurations uniformly sampled on the hemisphere surface.To compare our method with gradientfree methods, we also solve the task with Reinforcement Learning (RL), which has been widely used to train complex neural network policies for robotassisted dressing problems (Clegg et al., 2018). Specifically, we compare our gradientbased method with PPO (Schulman et al., 2017), a stateoftheart RL algorithm. We similarly sample a starting position from the hemisphere in each iteration when training PPO. For a fair comparison, we design the reward function to be the sum of the negative counterpart of plus a constant to avoid negative rewards, and observation of the environment to be . We evaluate Adam and PPO using as the common metric and plot the optimization curve at the bottom of Fig. 12. We see that both methods reach a similar final loss, but with our differentiable simulation framework, the gradientbased method reaches its final loss with a 85x speedup (Adam uses 23200 time steps; PPO uses 1978000 time steps).
After the training process converges, both our gradientbased method and PPO successfully move the hat onto the head from all 20 testing positions. We visualize the trajectories generated by Adam’s trained policy from 10 testing starting positions at the top of Fig. 12. Unlike existing differentiable simulation papers (Du et al., 2021; Hu et al., 2019, 2020; Liang et al., 2019) that train a closedloop network controller only for a fixed state, we highlight that we use differentiable simulation to train the network from multiple, random states and study its generalizability in a test set of unseen states. This allows us to conduct a fairer comparison between gradientbased optimization method and RL methods, which is typically overlooked in previous papers.
7. Conclusions, Limitation, and Future Work
In this paper, we presented a differentiable cloth simulator built on PD with SignoriniCoulomb frictional contact. Our differentiable simulator is different from existing papers in its simultaneous accommodation of rich and frequent (self)contact, SignoriniCoulomb contact law, and differentiability in cloth simulation. We analyzed the numerical properties of gradients from our differentiable simulator, including the source of discontinuities and its empirical speedup over gradientfree approaches in highdimensional problems. We additionally presented an iterative solver that exploits the contact gradients to speed up the backpropagation and observed a substantial speedup (up to an order of magnitude with lowprecision simulation) over direct solvers. Our differentiable cloth simulator enabled gradientbased optimization methods in a diverse set of applications, for which traditional gradientfree methods are generally much less sampling efficient. In particular, we presented a preliminary study on training a generalizable closedloop controller using differentiable simulation, in which our approach and PPO achieved comparable performance and generalizability, but we used much short amount of time.
There are still quite a few limitations in our method that are worth further investigation. First, since our method is built on PD, it also limits the choice of material models. It would be useful to generalize the current framework and support more physically accurate cloth material models, e.g., the piecewise linear elastic model described in Wang et al. (2011).
Second, the contact model we use from Ly et al. (2020) does not take into account vertexface or edgeedge selfcollisions. While we empirically observed that handling only vertexvertex selfcollisions managed to produce plausible results with medium mesh resolution, vertexface and edgeedge collision detection and handling is still highly desirable for a more physically realistic cloth simulator.
Third, although we have identified some possible sources of discontinuities and nonsmoothness in our differentiable simulator with empirical experiments, their effects on gradientbased optimization still require a thorough investigation. In particular, the locally bumpy energy landscape we observed in Fig. 3 due to changes in contact sets makes us question both the necessity and the usefulness of exact gradients in optimization, although Sec. 6 implies that gradients were still helpful in many downstream applications. Noting that the bumpiness in Fig. 3 is local and the global view of the energy landscape is still smooth, we hypothesize an inexact but smoothed gradients would be more powerful, which we leave as future work.
Fourth, there is no theoretical guarantee on the convergence of the iterative solver we implemented in backpropagation. Although nonconvergence is uncommon in our experiments, it introduces a costly switch to the slower direct solver, which we hope to fully resolve in future work.
Finally, many of our applications were in simulation only. It would be more exciting if these results could be replicated in realworld settings. We consider connecting our differentiable cloth simulator to more realworld applications, including realworld robotassisted dressing, material parameter identification for realworld fabric samples, computational design for sports suits, and so on. Closing the simtoreal gap is a nontrivial problem, in which we believe our differentiable simulator could play a beneficial role.
Appendix A Experiment Details
We provide detailed information for the examples shown in Sec. 6, including their setup, the exact form of their loss function, and their decision variables in optimization.
a.1. System Identification
Tshirt
The loss function is defined as
(31) 
where is the total number of time steps, the position of the mesh during optimization, and the position of the groundtruth simulation generated using our system. There are 4 parameters (6 DoFs) to be optimized: , where controls the stretching stiffness of the cloth material and controls the parameterized external wind force: each node receives a threedimensional wind force .
Sphere
The loss function is defined the same as in the “Tshirt” example above except that the total time . The decision variable is the 1DoF frictional coefficient of the sphere.
a.2. RobotAssisted Dressing
Hat
The loss function is defined as the L2 norm between the lasttimestep position of the hat and a target position generated by translating the initial pose of the hat onto the head:
(32) 
where is the total number of time steps. The trajectory of each of the two end effectors is controlled by a cubic Hermite spline. Each spline has parameters (9 DoFs): , where and are the two tangents of the spline and is the endpoint of the spline.
Sock
The goal is to optimize for the trajectory of the three end effectors holding onto a sock so that the sock can be put onto a synthetic foot model from an initial starting position. To guide the end effectors to first hook the opening of the sock onto the tip of the foot, then slide the sock upward onto the leg, the loss function is designed as
(33) 
where defines a set of manually selected keypoint correspondence pairs between the sock mesh and the foot model at halfway and the end of the simulation , respectively. We sum up the L2 norm of the position difference between each correspondence. The optimization parameters are the Hermite spline parameters for each of the four end effectors defined similarly as in the ”Hat” example (36 DoFs in total).
a.3. Inverse Design
Dress
The loss function is defined as
(34) 
where is the calculated target height of the bottom of the dress when the desired apex angle is reached. is the set of all points located at the bottom of the dress, and is the height of the point (corresponding to the y coordinate of the point in our implementation). There are 2 parameters (2 DoFs) that are being optimized: where is the density of the fabric and is the bending stiffness of the fabric.
a.4. RealtoSim Example
In this task, we optimize for the material parameters of the flag and the parameters of a simple wind model to match the motion trajectory of a flag to realworld captured data. The wind model is defined as
(35) 
where is a 3D vector to be optimized and a constant coefficient vector with each entry equal to a node’s inverse distance to the flag center evaluated at the first time step. There are 8 parameters to be optimized , where and are the stretching and bending stiffness of the fabric, the density of the fabric, and the parameters of the wind model defined above.
a.5. Hat Controller
In this task, the goal is to optimize for the neural network parameters of the hat controller so that the two end effectors put a hat on a head model from any starting position defined on a fixedradius hemisphere centered at the head model. During training, we uniformly sample 20 starting positions on the hemisphere for each epoch, and compute the loss averaged from all simulation sequences. We define the loss function as
(36) 
minimizes the distance change between the two end effectors and is defined as where and are the indices of the nodes pulled by the two end effectors. minimizes the L2distance between the poses of the hat and the target pose at the last few frames (20 in our implementation) and is defined as
(37) 
minimizes the orientation difference between the lasttimestep pose and the target pose of the hat and is defined as
(38) 
where is the up vector of the hat at the target pose and is defined similarly to the hat at the last time step.
References
 (1)
 Baraff and Witkin (1998) David Baraff and Andrew Witkin. 1998. Large Steps in Cloth Simulation. In Proceedings of the 25th Annual Conference on Computer Graphics and Interactive Techniques (SIGGRAPH ’98). Association for Computing Machinery, New York, NY, USA, 43–54.
 Bartle et al. (2016) Aric Bartle, Alla Sheffer, Vladimir G. Kim, Danny M. Kaufman, Nicholas Vining, and Floraine Berthouzoz. 2016. PhysicsDriven Pattern Adjustment for Direct 3D Garment Editing. ACM Trans. Graph. 35, 4, Article 50 (July 2016), 11 pages.
 Bouaziz et al. (2014) Sofien Bouaziz, Sebastian Martin, Tiantian Liu, Ladislav Kavan, and Mark Pauly. 2014. Projective Dynamics: Fusing Constraint Projections for Fast Simulation. ACM Trans. Graph. 33, 4, Article 154 (July 2014), 11 pages. https://doi.org/10.1145/2601097.2601116
 Bridson et al. (2002) Robert Bridson, Ronald Fedkiw, and John Anderson. 2002. Robust Treatment of Collisions, Contact and Friction for Cloth Animation. ACM Trans. Graph. 21, 3 (July 2002), 594–603.
 Bridson et al. (2003) R. Bridson, S. Marino, and R. Fedkiw. 2003. Simulation of Clothing with Folds and Wrinkles. In Proceedings of the 2003 ACM SIGGRAPH/Eurographics Symposium on Computer Animation (San Diego, California) (SCA ’03). Eurographics Association, Goslar, DEU, 28–36.
 Brogliato (2016) Bernard Brogliato. 2016. Nonsmooth Lagrangian Systems. Springer International Publishing, Cham, 241–370.
 Brouet et al. (2012) Remi Brouet, Alla Sheffer, Laurence Boissieux, and MariePaule Cani. 2012. Design Preserving Garment Transfer. ACM Trans. Graph. 31, 4, Article 36 (July 2012), 11 pages.
 Brown et al. (2018) George E. Brown, Matthew Overby, Zahra Forootaninia, and Rahul Narain. 2018. Accurate Dissipative Forces in Optimization Integrators. ACM Trans. Graph. 37, 6, Article 282 (Dec. 2018), 14 pages.
 Choi and Ko (2002) KwangJin Choi and HyeongSeok Ko. 2002. Stable but Responsive Cloth. ACM Trans. Graph. 21, 3 (July 2002), 604–611.
 Clegg et al. (2020) Alexander Clegg, Zackory Erickson, Patrick Grady, Greg Turk, Charles C. Kemp, and C. Karen Liu. 2020. Learning to Collaborate From Simulation for RobotAssisted Dressing. IEEE Robotics and Automation Letters 5, 2 (2020), 2746–2753.
 Clegg et al. (2018) Alexander Clegg, Wenhao Yu, Jie Tan, C. Karen Liu, and Greg Turk. 2018. Learning to Dress: Synthesizing Human Dressing Motion via Deep Reinforcement Learning. ACM Trans. Graph. 37, 6, Article 179 (Dec. 2018), 10 pages.
 Clyde et al. (2017) David Clyde, Joseph Teran, and Rasmus Tamstorf. 2017. Modeling and DataDriven Parameter Estimation for Woven Fabrics. In Proceedings of the ACM SIGGRAPH / Eurographics Symposium on Computer Animation (Los Angeles, California) (SCA ’17). Association for Computing Machinery, New York, NY, USA, Article 17, 11 pages.
 Dario Bellicoso et al. (2016) C. Dario Bellicoso, Christian Gehring, Jemin Hwangbo, Péter Fankhauser, and Marco Hutter. 2016. Perceptionless terrain adaptation through whole body control and hierarchical optimization. In 2016 IEEERAS 16th International Conference on Humanoid Robots (Humanoids). 558–564. https://doi.org/10.1109/HUMANOIDS.2016.7803330
 de Avila BelbutePeres et al. (2018a) Filipe de Avila BelbutePeres, Kevin Smith, Kelsey Allen, Josh Tenenbaum, and J Zico Kolter. 2018a. Endtoend differentiable physics for learning and control. Advances in neural information processing systems 31 (2018), 7178–7189.
 de Avila BelbutePeres et al. (2018b) Filipe de Avila BelbutePeres, Kevin Smith, Kelsey Allen, Josh Tenenbaum, and Zico Kolter. 2018b. EndtoEnd Differentiable Physics for Learning and Control. In Advances in Neural Information Processing Systems. 7178–7189.
 Degrave et al. (2019) Jonas Degrave, Michiel Hermans, Joni Dambre, and Francis wyffels. 2019. A Differentiable Physics Engine for Deep Learning in Robotics. Frontiers in Neurorobotics 13 (2019), 6.
 Du et al. (2021) Tao Du, Kui Wu, Pingchuan Ma, Sebastien Wah, Andrew Spielberg, Daniela Rus, and Wojciech Matusik. 2021. DiffPD: Differentiable Projective Dynamics with Contact. arXiv:2101.05917 [cs.LG]
 Du et al. (2020) Tao Du, Kui Wu, Andrew Spielberg, Wojciech Matusik, Bo Zhu, and Eftychios Sifakis. 2020. Functional Optimization of Fluidic Devices with Differentiable Stokes Flow. ACM Trans. Graph. 39, 6, Article 197 (Dec. 2020), 15 pages. https://doi.org/10.1145/3414685.3417795
 Fratarcangeli et al. (2016) Marco Fratarcangeli, Valentina Tibaldo, and Fabio Pellacini. 2016. Vivace: A Practical GaussSeidel Method for Stable Soft Body Dynamics. ACM Trans. Graph. 35, 6, Article 214 (Nov. 2016), 9 pages.
 Geilinger et al. (2020) Moritz Geilinger, David Hahn, Jonas Zehnder, Moritz Bächer, Bernhard Thomaszewski, and Stelian Coros. 2020. ADD: analytically differentiable dynamics for multibody systems with frictional contact. ACM Transactions on Graphics (TOG) 39, 6 (2020), 1–15.
 Guan et al. (2012) Peng Guan, Loretta Reiss, David A. Hirshberg, Alexander Weiss, and Michael J. Black. 2012. DRAPE: DRessing Any PErson. ACM Trans. Graph. 31, 4, Article 35 (July 2012), 10 pages.
 Guennebaud et al. (2010) Gaël Guennebaud, Benoît Jacob, et al. 2010. Eigen v3. http://eigen.tuxfamily.org.
 Hahn et al. (2019) David Hahn, Pol Banzet, James M Bern, and Stelian Coros. 2019. Real2sim: Viscoelastic parameter estimation from dynamic motion. ACM Transactions on Graphics (TOG) 38, 6 (2019), 1–13.
 Hansen (2006) Nikolaus Hansen. 2006. The CMA Evolution Strategy: A Comparing Review.
 Harmon et al. (2008) David Harmon, Etienne Vouga, Rasmus Tamstorf, and Eitan Grinspun. 2008. Robust Treatment of Simultaneous Collisions. In ACM SIGGRAPH 2008 Papers (Los Angeles, California) (SIGGRAPH ’08). Association for Computing Machinery, New York, NY, USA, Article 23, 4 pages.
 Holl et al. (2020) Philipp Holl, Nils Thuerey, and Vladlen Koltun. 2020. Learning to Control PDEs with Differentiable Physics. In International Conference on Learning Representations.
 Hu et al. (2020) Yuanming Hu, Luke Anderson, TzuMao Li, Qi Sun, Nathan Carr, Jonathan RaganKelley, and Frédo Durand. 2020. DiffTaichi: Differentiable Programming for Physical Simulation. International Conference on Learning Representations (2020).
 Hu et al. (2019) Yuanming Hu, Jiancheng Liu, Andrew Spielberg, Joshua B Tenenbaum, William T Freeman, Jiajun Wu, Daniela Rus, and Wojciech Matusik. 2019. Chainqueen: A realtime differentiable physical simulator for soft robotics. In 2019 International conference on robotics and automation (ICRA). IEEE, 6265–6271.
 Kang et al. (2021) Dongho Kang, Simon Zimmermann, and Stelian Coros. 2021. Animal Gaits on Quadrupedal Robots Using Motion Matching and ModelBased Control. In 2021 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS). IEEE.
 Kingma and Ba (2017) Diederik P. Kingma and Jimmy Ba. 2017. Adam: A Method for Stochastic Optimization. arXiv:1412.6980 [cs.LG]
 Komaritzan and Botsch (2018) Martin Komaritzan and Mario Botsch. 2018. Projective Skinning. Proc. ACM Comput. Graph. Interact. Tech. 1, 1, Article 12 (July 2018), 19 pages. https://doi.org/10.1145/3203203
 Li et al. (2020b) Cheng Li, Min Tang, Ruofeng Tong, Ming Cai, Jieyi Zhao, and Dinesh Manocha. 2020b. PCloth: Interactive Complex Cloth Simulation on MultiGPU Systems Using Dynamic Matrix Assembly and Pipelined Implicit Integrators. ACM Trans. Graph. 39, 6, Article 180 (Nov. 2020), 15 pages.
 Li et al. (2018b) Jie Li, Gilles Daviet, Rahul Narain, Florence BertailsDescoubes, Matthew Overby, George E. Brown, and Laurence Boissieux. 2018b. An Implicit Frictional Contact Solver for Adaptive Cloth Simulation. ACM Trans. Graph. 37, 4, Article 52 (July 2018), 15 pages.
 Li et al. (2020a) Minchen Li, Zachary Ferguson, Teseo Schneider, Timothy Langlois, Denis Zorin, Daniele Panozzo, Chenfanfu Jiang, and Danny M. Kaufman. 2020a. Incremental Potential Contact: Intersectionand InversionFree, LargeDeformation Dynamics. ACM Trans. Graph. 39, 4, Article 49 (July 2020), 20 pages.
 Li et al. (2021) Minchen Li, Danny M. Kaufman, and Chenfanfu Jiang. 2021. Codimensional Incremental Potential Contact. ACM Trans. Graph. (SIGGRAPH) 40, 4, Article 170 (2021).
 Li et al. (2018c) Minchen Li, Alla Sheffer, Eitan Grinspun, and Nicholas Vining. 2018c. Foldsketch: Enriching Garments with Physically Reproducible Folds. ACM Trans. Graph. 37, 4, Article 133 (July 2018), 13 pages.
 Li et al. (2018a) TzuMao Li, Miika Aittala, Frédo Durand, and Jaakko Lehtinen. 2018a. Differentiable Monte Carlo Ray Tracing through Edge Sampling. ACM Trans. Graph. (Proc. SIGGRAPH Asia) 37, 6 (2018), 222:1–222:11.
 Li et al. (2019) Yunzhu Li, Jiajun Wu, Russ Tedrake, Joshua Tenenbaum, and Antonio Torralba. 2019. Learning Particle Dynamics for Manipulating Rigid Bodies, Deformable Objects, and Fluids. In International Conference on Learning Representations.
 Liang et al. (2019) Junbang Liang, Ming Lin, and Vladlen Koltun. 2019. Differentiable Cloth Simulation for Inverse Problems. In Advances in Neural Information Processing Systems, H. Wallach, H. Larochelle, A. Beygelzimer, F. d'AlchéBuc, E. Fox, and R. Garnett (Eds.), Vol. 32. Curran Associates, Inc. https://proceedings.neurips.cc/paper/2019/file/28f0b864598a1291557bed248a998d4ePaper.pdf
 Liu and Nocedal (1989) Dong C. Liu and Jorge Nocedal. 1989. On the Limited Memory BFGS Method for Large Scale Optimization. MATHEMATICAL PROGRAMMING 45 (1989), 503–528.
 Liu et al. (2013) Tiantian Liu, Adam W. Bargteil, James F. O’Brien, and Ladislav Kavan. 2013. Fast Simulation of MassSpring Systems. ACM Trans. Graph. 32, 6, Article 214 (Nov. 2013), 7 pages.
 Liu et al. (2017) Tiantian Liu, Sofien Bouaziz, and Ladislav Kavan. 2017. Quasinewton methods for realtime simulation of hyperelastic materials. Acm Transactions on Graphics (TOG) 36, 3 (2017), 1–16.
 Ly et al. (2018) Mickaël Ly, Romain Casati, Florence BertailsDescoubes, Mélina Skouras, and Laurence Boissieux. 2018. Inverse Elastic Shell Design with Contact and Friction. ACM Trans. Graph. 37, 6, Article 201 (Dec. 2018), 16 pages.
 Ly et al. (2020) Mickaël Ly, Jean Jouve, Laurence Boissieux, and Florence BertailsDescoubes. 2020. Projective Dynamics with Dry Frictional Contact. ACM Trans. Graph. 39, 4, Article 57 (July 2020), 8 pages.
 Macklin et al. (2020) M. Macklin, K. Erleben, M. Müller, N. Chentanez, S. Jeschke, and T. Y. Kim. 2020. Primal/Dual Descent Methods for Dynamics. In Proceedings of the ACM SIGGRAPH/Eurographics Symposium on Computer Animation (Virtual Event, Canada) (SCA ’20). Eurographics Association, Goslar, DEU, Article 9, 12 pages.
 Macklin et al. (2016) Miles Macklin, Matthias Müller, and Nuttapong Chentanez. 2016. XPBD: PositionBased Simulation of Compliant Constrained Dynamics. In Proceedings of the 9th International Conference on Motion in Games (Burlingame, California) (MIG ’16). Association for Computing Machinery, New York, NY, USA, 49–54.
 Martin et al. (2011) Sebastian Martin, Bernhard Thomaszewski, Eitan Grinspun, and Markus Gross. 2011. ExampleBased Elastic Materials. ACM Trans. Graph. 30, 4, Article 72 (July 2011), 8 pages. https://doi.org/10.1145/2010324.1964967
 McNamara et al. (2004) Antoine McNamara, Adrien Treuille, Zoran Popović, and Jos Stam. 2004. Fluid Control Using the Adjoint Method. ACM Trans. Graph. 23, 3 (Aug. 2004), 449–456.
 Miguel et al. (2012) E. Miguel, D. Bradley, B. Thomaszewski, B. Bickel, W. Matusik, M. A. Otaduy, and S. Marschner. 2012. DataDriven Estimation of Cloth Simulation Models. Comput. Graph. Forum 31, 2pt2 (May 2012), 519–528.
 Miguel et al. (2013) Eder Miguel, Rasmus Tamstorf, Derek Bradley, Sara C. Schvartzman, Bernhard Thomaszewski, Bernd Bickel, Wojciech Matusik, Steve Marschner, and Miguel A. Otaduy. 2013. Modeling and Estimation of Internal Friction in Cloth. ACM Trans. Graph. 32, 6, Article 212 (Nov. 2013), 10 pages.
 Mistry et al. (2010) Michael Mistry, Jonas Buchli, and Stefan Schaal. 2010. Inverse dynamics control of floating base systems using orthogonal decomposition. In 2010 IEEE International Conference on Robotics and Automation. 3406–3412. https://doi.org/10.1109/ROBOT.2010.5509646
 Montes et al. (2020) Juan Montes, Bernhard Thomaszewski, Sudhir Mudur, and Tiberiu Popa. 2020. Computational Design of Skintight Clothing. ACM Trans. Graph. 39, 4, Article 105 (July 2020), 12 pages.
 Moré and Thuente (1994) Jorge J. Moré and David J. Thuente. 1994. Line Search Algorithms with Guaranteed Sufficient Decrease. ACM Trans. Math. Softw. 20, 3 (Sept. 1994), 286–307. https://doi.org/10.1145/192115.192132
 Müller et al. (2007) Matthias Müller, Bruno Heidelberger, Marcus Hennix, and John Ratcliff. 2007. Position Based Dynamics. J. Vis. Comun. Image Represent. 18, 2 (April 2007), 109–118.
 Murthy et al. (2021) J. Krishna Murthy, Miles Macklin, Florian Golemo, Vikram Voleti, Linda Petrini, Martin Weiss, Breandan Considine, Jérôme ParentLévesque, Kevin Xie, Kenny Erleben, Liam Paull, Florian Shkurti, Derek Nowrouzezahrai, and Sanja Fidler. 2021. gradSim: Differentiable simulation for system identification and visuomotor control. In International Conference on Learning Representations. https://openreview.net/forum?id=c_E8kFWfhp0
 Narain et al. (2012) Rahul Narain, Armin Samii, and James F. O’Brien. 2012. Adaptive Anisotropic Remeshing for Cloth Simulation. ACM Trans. Graph. 31, 6, Article 152 (Nov. 2012), 10 pages.
 Otaduy et al. (2009) Miguel A. Otaduy, Rasmus Tamstorf, Denis Steinemann, and Markus Gross. 2009. Implicit Contact Handling for Deformable Objects. Computer Graphics Forum 28, 2 (2009), 559–568.
 Overby et al. (2017) Matthew Overby, George E. Brown, Jie Li, and Rahul Narain. 2017. ADMM Projective Dynamics: Fast Simulation of Hyperelastic Models with Dynamic Constraints. IEEE Transactions on Visualization and Computer Graphics 23, 10 (Oct 2017), 2222–2234.
 Popović et al. (2003) Jovan Popović, Steven Seitz, and Michael Erdmann. 2003. Motion Sketching for Control of RigidBody Simulations. ACM Trans. Graph. 22, 4 (Oct. 2003), 1034–1054.
 Popović et al. (2000) Jovan Popović, Steven M. Seitz, Michael Erdmann, Zoran Popović, and Andrew Witkin. 2000. Interactive Manipulation of Rigid Body Simulations. In Proceedings of the 27th Annual Conference on Computer Graphics and Interactive Techniques (SIGGRAPH ’00). ACM Press/AddisonWesley Publishing Co., USA, 209–217. https://doi.org/10.1145/344779.344880
 Provot (1997) Xavier Provot. 1997. Collision and selfcollision handling in cloth model dedicated to design garments. In Computer Animation and Simulation ’97, Daniel Thalmann and Michiel van de Panne (Eds.). Springer Vienna, Vienna, 177–189.
 Qiao et al. (2020) YiLing Qiao, Junbang Liang, Vladlen Koltun, and Ming Lin. 2020. Scalable Differentiable Physics for Learning and Control. In International Conference on Machine Learning.
 SanchezGonzalez et al. (2020) Alvaro SanchezGonzalez, Jonathan Godwin, Tobias Pfaff, Rex Ying, Jure Leskovec, and Peter Battaglia. 2020. Learning to Simulate Complex Physics with Graph Networks. In International Conference on Machine Learning.
 Schenck and Fox (2018) Connor Schenck and Dieter Fox. 2018. SPNets: Differentiable Fluid Dynamics for Deep Neural Networks. Conference on Robot Learning (CoRL) (2018).
 Schulman et al. (2017) John Schulman, Filip Wolski, Prafulla Dhariwal, Alec Radford, and Oleg Klimov. 2017. Proximal policy optimization algorithms. arXiv preprint arXiv:1707.06347 (2017).
 Stern and Desbrun (2006) Ari Stern and Mathieu Desbrun. 2006. Discrete Geometric Mechanics for Variational Time Integrators. In ACM SIGGRAPH 2006 Courses (Boston, Massachusetts) (SIGGRAPH ’06). Association for Computing Machinery, New York, NY, USA, 75–80.
 Stern and Grinspun (2009) Ari Stern and Eitan Grinspun. 2009. ImplicitExplicit Variational Integration of Highly Oscillatory Problems. Multiscale Modeling & Simulation 7, 4 (2009), 1779–1794.
 Tamstorf et al. (2015) Rasmus Tamstorf, Toby Jones, and Stephen F. McCormick. 2015. Smoothed Aggregation Multigrid for Cloth Simulation. ACM Trans. Graph. 34, 6, Article 245 (Oct. 2015), 13 pages.
 Terzopoulos et al. (1987) Demetri Terzopoulos, John Platt, Alan Barr, and Kurt Fleischer. 1987. Elastically Deformable Models. SIGGRAPH Comput. Graph. 21, 4 (Aug. 1987), 205–214.
 Toussaint et al. (2018) Marc Toussaint, Kelsey Allen, Kevin Smith, and Joshua Tenenbaum. 2018. Differentiable Physics and Stable Modes for ToolUse and Manipulation Planning. In Robotics: Science and Systems, Vol. 2.
 Treuille et al. (2003) Adrien Treuille, Antoine McNamara, Zoran Popović, and Jos Stam. 2003. Keyframe Control of Smoke Simulations. ACM Trans. Graph. 22, 3 (July 2003), 716–723.
 Umetani et al. (2011) Nobuyuki Umetani, Danny M. Kaufman, Takeo Igarashi, and Eitan Grinspun. 2011. Sensitive Couture for Interactive Garment Modeling and Editing. ACM Trans. Graph. 30, 4, Article 90 (July 2011), 12 pages.
 Wang (2015) Huamin Wang. 2015. A Chebyshev SemiIterative Approach for Accelerating Projective and PositionBased Dynamics. ACM Trans. Graph. 34, 6, Article 246 (Oct. 2015), 9 pages.
 Wang (2018) Huamin Wang. 2018. RuleFree Sewing Pattern Adjustment with Precision and Efficiency. ACM Trans. Graph. 37, 4, Article 53 (July 2018), 13 pages.
 Wang et al. (2011) Huamin Wang, James F. O’Brien, and Ravi Ramamoorthi. 2011. DataDriven Elastic Models for Cloth: Modeling and Measurement. In ACM SIGGRAPH 2011 Papers (Vancouver, British Columbia, Canada) (SIGGRAPH ’11). Association for Computing Machinery, New York, NY, USA, Article 71, 12 pages.
 Wang and Yang (2016) Huamin Wang and Yin Yang. 2016. Descent Methods for Elastic Body Simulation on the GPU. ACM Trans. Graph. 35, 6, Article 212 (Nov. 2016), 10 pages.
 Wang et al. (2018) Zhendong Wang, Longhua Wu, Marco Fratarcangeli, Min Tang, and Huamin Wang. 2018. Parallel Multigrid for Nonlinear Cloth Simulation. Computer Graphics Forum 37, 7 (2018), 131–141.
 Werling et al. (2021) Keenon Werling, Dalton Omens, Jeongseok Lee, Ioannis Exarchos, and C. Karen Liu. 2021. Fast and FeatureComplete Differentiable Physics Engine for Articulated Rigid Bodies with Contact Constraints. In Proceedings of Robotics: Science and Systems. Virtual. https://doi.org/10.15607/RSS.2021.XVII.034
 White et al. (2007) Ryan White, Keenan Crane, and David Forsyth. 2007. Capturing and Animating Occluded Cloth. In ACM Transactions on Graphics (SIGGRAPH).
 Witkin and Kass (1988) Andrew Witkin and Michael Kass. 1988. Spacetime Constraints. In Proceedings of the 15th Annual Conference on Computer Graphics and Interactive Techniques (SIGGRAPH ’88). Association for Computing Machinery, New York, NY, USA, 159–168. https://doi.org/10.1145/54852.378507
 Wojtan et al. (2006) Chris Wojtan, Peter Mucha, and Greg Turk. 2006. Keyframe Control of Complex Particle Systems Using the Adjoint Method. In Proceedings of the 2006 ACM SIGGRAPH/Eurographics Symposium on Computer Animation (Vienna, Austria) (SCA ’06). Eurographics Association, Goslar, DEU, 15–23.
 Xian et al. (2019) Zangyueyang Xian, Xin Tong, and Tiantian Liu. 2019. A Scalable Galerkin Multigrid Method for RealTime Simulation of Deformable Objects. ACM Trans. Graph. 38, 6, Article 162 (Nov. 2019), 13 pages.
 Xu et al. (2021) Jie Xu, Tao Chen, Lara Zlokapa, Michael Foshey, Wojciech Matusik, Shinjiro Sueda, and Pulkit Agrawal. 2021. An EndtoEnd Differentiable Framework for ContactAware Robot Design. In Proceedings of Robotics: Science and Systems. Virtual. https://doi.org/10.15607/RSS.2021.XVII.008
 Yixuan (2021) Qiu Yixuan. 2021. LBFGS++. https://github.com/yixuan/LBFGSpp/.
Comments
There are no comments yet.