CS384G Final Project: Realistic Cloth Simulation

December 9, 2004


Jenn Sartor
Peter Djeu

splash image

I.  Overview

Realistic cloth simulation is a computational expensive procedure.  We decided to pursue this topic with the hope of understanding the problem domain while finding new ways to apply many of the concepts we have learned throughout this semester.

We simulate cloth in an animation environment using a particle system and a two-dimensional grid of particles, which is visualized as a triangle mesh.  We compute per-vertex normals and apply a rug pattern texture map for added realism.  Our simulator is an extension of the Project 3 Animator assignment.


II. Simulating Cloth Dynamics

II. Part A.  Forces

The most important component of simulating cloth is accurately computing the internal and external forces acting upon the cloth particles.  The forces we cover include:

We will explain each of the forces in more detail below.

type 1 type 2

Type I
Type I forces model the stretching and sheering of cloth.  This is a spring force between a particle and its 8 nearest neighbors.  The red edges in Figure 2 (from [Choi 2002]) depict all type I interactions for the center particle.  In other words, each particle is connect by a type I force to all particles that are within one hop in the grid.  In general this allows the cloth to stretch and sheer, but only to a limited degree.  Type I interaction is governed by two constants: 1) rest length (L), which specified the distance past which a spring will exert contraction force and 2) the spring constant (k_s), which is directly proportional to the amount of recovery force the spring exerts when it is stetched past its rest length.

The outer edge of the cloth presents a problem due to the assymetric spring forces on the particles here.  In particular, each of these particles has either 3 or 5 neighbors (rather than the balanced grid of 8 neighbors), so the type I force exhibited by the interior neighbors will uniformly pull the outermost particles towards the center of the cloth.  In general, this assymetry causes the cloth to uniformly shrink towards its center.  To compensate for this problem, we increase the mass of these outer particles tenfold.

Another property of type I is that it only models the cloth's resistance to stretching, which means that if the distance between two particles i and j with type I interaction is less than L, no type I force is exerted between these particles.  This is different from a traditional spring under Hooke's law, where a recovery force would be exerted if two particles were less than distance L apart.


Type II
Type II forces model the repulsion force when cloth buckles (i.e. when cloth collapses in on itself).  These interactions are modeled between a particle and the 8 neighbors that are two particles away.  The type II interactions are shown in blue in the figure above.  In the 1-D case, this force creates realistic cloth buckling by maintaining a certain angle theta between the two grid edges that leave a particle.  In Figure 4 (from [Choi 2002]), the particle endpoints of the line segment represent the two particles on which the type II interaction is taking affect.  The particle that they span is the point at which theta is being regulated.

type2_details


Type II force is calculated using a flexural rigidity constant, k_b, that can vary the strength of repulsion.  This calculation also uses the inverse sinc function ([Choi 2002]) which we approximate with a finite series. The type II force is not exerted between particles if the distance between the two particles in the interaction is greater than or equal to the spring rest length L.


Damping
We use a linear damping force between all pairs (particle i, particle j) which have either a type I or type II interaction.  The damping force is a constant factor of  the difference of their velocities.  This force stops infinite in-plane oscillations normally caused by spring forces.   The strength of this force is changed with the k_d damping constant.


Gravity
We modeled gravity as a constant force across each time step (9.8 m/s2 downwards).


Wind
For wind simulations, we modeled a wind force that only acts in the positive and negative x directions.  In order to model variable wind accurately, we multiply the speed by the sin function evaluated with the current time in the particle simulation.  In order for the wind to vary at different particles (which prevents all particles from moving as a wall), we vary the speed based on a particle's distance from the center of the grid.  We found the most appealing affect involved making the wind force strongest at the center of the cloth.  The wind strength is also varied by the angle at which it hits the cloth (strongest when the wind force coincides with the cloth normal).  To anchor the cloth, three particle positions are held constant along the top row, which models the cloth pinned to a clothesline.

II. Part B.  Numerical Solvers

Euler's Method
For all movies shown on this page, we used Euler's method along with force accumulators to update particle position and velocity, in much the same way as the particle system was updated in Project 3.  This method is well-known.  To avoid numerical instabilities, we capped the maximum allowable particle velocity and divided each time step into many smaller increments which we call mini time steps.  Euler's method becomes more unstable with larger time steps, so the mini time steps cut down on mathematical error.

Preconditioned Conjugate Gradient Method
In order to reduce numerical instabilities, we implemented the preconditioned conjugate gradient (PCG) method for updating particles.  This method is more stable than Euler's method, and instead of using force accumulators, it uses the particle positions and velocities over the last two time steps -- as well as the particle forces at the last time step -- to determine current particle positions, velocities, and forces.

Under the PCG method, each force described in Section I has a corresponding Jacobian matrix version that is used in place of the force vector.  We use Jacobian matrices for all of the forces listed above (with the exception of the wind force).
Each of these Jacobians is a three by three matrix and belongs within a much larger partial derivative matrix ðf/ðx.  Each of the these partial derivatives that vary with i and j are entered into the ðf/ðx matrix, which incidentally is a 3N by 3N matrix where N is the number of particles being simulated.  So the Jacobian of particle i with respect to particle j would be start at row i * 3 and column j * 3.

Using the partial derivative matrix in place of the force accumulators, as well as the particle masses, the time step, and the previous particle positions, velocities, and forces, we compute the system matrix A in the form Ax = b (see Equation (20) below [Choi 2002]).  A is a 3N by 3N matrix, and x and b are vectors of length 3N.   A and b are  known, and x is solved for using the PCG solver.  The x vector represents the positions  of all particles at the current time step (what we are solving for) minus the positions of all particles at the previous time step.  Therefore, x = xn - xn-1.   In Equation (20), the vector x (length 3N) that we are solving for appears at the end of the first line.  After we use PCG to solve this system for x, we can solve for the vector of updated particle positions ( xn ).  From there we can compute the vector of updated particle velocities and the vector of updated forces as well.  We used and slightly modified the code for PCG from Numerical Recipes, a numerical analysis programming book [Press 1996] .

Unfortunately the PCG method is computationally expensive and uses a lot of memory.  Even with very small simulations, we were not able to get the PCG method to model cloth accurately.
cloth_eqns

III. Collision Detection

We include two type of intersectable shapes in our cloth simulation: spheres and axis-aligned boxes.  Our scene can include any number of these intersectable shapes, and via collision detection, the cloth will react realistically when it is the vicinity of these shapes.  At each time step, every particle is intersected against every shape in the scene as well as the floor.

III. Part A.  Axis-Aligned Boxes

Collision detection for a box entailes performing collision detection with  the six clipped planes that comprise the box.  We perform particle-clipped plane collision detection by testing if the previous particle position was on one side of the plane, and the current particle position is on the other side of the plane (which is undesired).  When a collision is detected, the coordinates for two of the three dimensions (x and z) are reused from the previous time step, and the y coordinate is calculated so that the particle lies on the proper side of the clipped plane.  This algorithm was chosen because most collisions occur due to gravity, so the y dimension is the most important dimension to correct.

III. Part B.  Spheres

For particle-sphere collision detection, we first compute the distance between a particle and the center of the sphere (in vector form).  If this distance is smaller than the radius of the sphere,  the distance vector is normalized, then scaled by  the radius so that the particle is placed on the outside of the sphere.

III. Part C.  Shape Visualization

In order to prevent our scene objects from peeking through the cloth due to boundary conditions, all geometric shapes were drawn slightly smaller (by an epsilon) than their actual dimensions.  An alternative solution is to perform collision detection with a shape that is epsilon larger than the actual objects we are drawing, but this solution turns out to be hairier to implement without yielding appreciably better results.


IV. Cloth Constants

As mentioned above, we can modify many constants in this cloth simulation.  Varying each of the constants can radically change the look and feel of our simluated cloth.  After talking to the authors of the cloth paper (Choi and Ko), we determined values that work well for the spring, flexural rigidity, and damping constants based on the the simulation they were needed in.  Each of these can be varied via a slider in the animation window.

The number of particles is fixed at compile time, but the distance between particles and the spring rest length can also be changed in sliders in the animation window.  Also for Euler's method, the number of mini time steps used, and the velocity cap can be set per-frame at runtime.


V. Results

V.  Part A.  Performance

The following compute times were needed for a 30 second simulation of dropping cloth onto a box (see Box drop in Artifacts section).  The simulation ran with 30 fps and 10 mini time steps for every frame.

Grid Resolution
 Compute Time
15 x 15
 45 seconds
30 x 30
 134 seconds
45 x 45
 402 seconds


V. Part B.  Artifacts


We used our cloth simulation to generate the following artifacts (download by clicking on the images).  Results were generated using a 15 x 15 grid of particles unless otherwise noted.  We hope you enjoy watching these movies!


Simple 1-D simulation - Demonstrates type I forces in a simple 1-D simulation, where each particle location is represented with a sphere.  A 5 x 1 cloth grid is used.  Notice that the type I forces eventually reach an equilibrium with the force of gravity.

5particle_thumb


Box drop - Cloth is dropped onto a box.  The flair-out effect of the corners of the cloth provides a particularly nice visual effect.

box_thumb


Sphere drop - Cloth is dropped onto a sphere.  The cloth slips off the face of the sphere due to the effect of gravity.

sphere_thumb


Box-Sphere drop - Cloth is dropped onto a sphere that is floating above a box.  Demonstrates multiple collidable objects in our scene.

box_sphere_thumb


Wind tunnel - Cloth is subjected to the wind, followed by a cannonball attack.

wind_thumb


Acknowledgements

We would like to kindly thank Kwang-Jin Choi and Hyeong-Seok Ko for sharing their patience and expertise during the development phase of our project.


Bibliography


Bridson, R.  Fedkiw, R.  and Anderson, J. 2002.  Robust treatment of collisions, contact and friction for cloth animation.    In Proceedings of the 29th Annual Conference on Computer Graphics and Interactive Techniques, 594 - 603.
link to paper


Choi, K.-J. and Ko, H.-S. 2002.  Stable but responsive cloth.    In Proceedings of the 29th Annual Conference on Computer Graphics and Interactive Techniques, 604--611.
link to paper


Press, W.  Teukolsky, S.  Vetterling, W.  Flannery, B. 1996 Reprint.  Numerical Recipes in C: The Art of Scientific Computing.  Cambridge University Press.