COMP 559 - Winter 2011 - Assignment 2
Rigid Body Collision and Contact

Due noon Friday February 11

Getting Started

This assignment is based on a classic Cornell assignment. It involves the simulation of contact in a system of 2D rigid bodies. There are three main parts: rotational dynamics, hierarchical collision detection, and a velocity level resolution of contact using a projected Gauss-Seidel complementarity constraint solver. The starter code provides you with a working simulation, except that there is no rotational dynamics, collision detection is done using a naive all-pairs approach, and frictionless contacts are simulated with penalty springs.

Download the provided code from WebCT, and dump it into a new java project. Do not change the package names as this will interfere with marking. You can create any new classes you like and modify the provided ones, but you'll likely want to focus your efforts in two areas, the RigidBody class and the Collision processor class. Read through all the provided code and javadoc to familiarize yourself with the structure. Here is a rough description of the different classes provided:

A2App - Main entry point for the application. You do not need to modify this code. It creates a swing interface for adjusting viewing and simulation parameters (explore these parameters). Note also 2 buttons to let you set the drawing canvas to an optimal size for recording frames to make a video. The keyboard interface provides the following:
- space - toggles running of the simulation
- s - steps the simulation
- r - reset
- c - clear
- l - load an image file with a file selector
- j - jiggle all bodies a small amount
- f - sets up the default factory for stress testing
- g - sets up a factory from an image loaded with the file selector
- left arrow - load previous image file in data directory
- right arrow - load next image file in data directory
- up arrow - increase substeps and adjust step size
- down arrow - decreaes substeps and adjust step size
- enter - toggles recording of canvas to stills directory
- esc - quit
Block - The Block class represents a single pixel in the orignal image. It also keeps track of the block position in the body frame.
BVNode - Bounding volume node used to build bounding volume trees.
CollisionProcessor - Class for detecting and resolving collisions. You will need to add heirarchical collision detection and a solver for contact constraints.
Contact - Implementation of a contact between two bodies. You will want to add code to compute the constraint Jacobian.
Disc - Circular disc for collision processing.
Factory - RigidBody factory for stress testing.
ImageBlocker - Creates rigid bodies from an image.
MinimumEnclosingCircle - Algorithm for computing Minimum Enclosing Circles using Welzl's algorithm.
MouseSpringForce - Mouse spring for interacting with rigid bodies.
RigidBody - Simple 2D rigid body based on image samples. You will need to add code to do rotational dynamics.
RigidBodySystem - Maintains a list of RigidBody objects, and provides methods for collision processing and numerical integration.
RigidTransform - A class for 2D rigid transformations.

The provided code zip file also includes a number of test images. Any new files you add should be in png format. Objects that are defined by connected components of non-white pixels (note there is a threshold so that near white pixels are also considered to be empty). Objects that consist entirely of shades of blue are pinned. Feel free to create your own files to create simple tests, or something spectacular. Share your images on WebCT for others to try!

The code makes use of a new mintools.jar, also included in the zip file. It has been updated for this assignemnt with a file selector class and a fix for the slider initialization problem. The source for the tools is provided for your convenience, but you should not need to make changes. JOGL for OpenGL, and vecmath are also used by the provided code, and you can reuse the libraries you downloaded for assignment 1. You can also use MTJ in this assignment if you choose, however, your code will probably be much faster if you don't allocate and assemble the sparse matrices needed in your Gauss Seidel solver.

Steps and Objectives (10/10)

The important parts of the program that you will need to modify or write to complete this assignment are labeled with TODO comments. When you are running large test cases you will want to use the JVM's -server option to encourage aggressive runtime optimizations. You may also find profiling useful to help identify and remove bottlenecks.

Rotational Rigid Body Dynamics (2 marks)

Most of the rigid body dynamics code is provided already, but only linear forces and dynamics are implemented. Modify the RigidBody class to compute torques caused by forces acting at a specified point, and use these torques to advance the state of the system (angular velocity and angular position).
Hierarchical Collision Detection (2 marks)

Collision detection is extremely slow because all blocks are checked with all others. Code to build a bounding disc tree is already provided for you, but it isn't used. Use the BVNode trees of pairs or rigid bodies to speed up body-body tests.
LCP Solve (5 marks)

Replace the penalty based contact response with a velocity level update to satisfy the contact constraints. You must implement a projected Gauss-Seidel solver as discussed in class [Erleben 2007]. As contact detection is already done for the penalty springs, you will already have a list of Contact objects where the normal is the direction between two overlapping blocks and the contact point is the midpoint between them. There are three parts to this objective.
1. Compute the contact Jacobian.
2. Set up your w = A lambda + b system of equations by computing the b vector, and write code for evaluation of Ai lambda (the ith row of A multiplied by lambda). To be efficient, you will want to do this without assembling A, and instead, you should implement a sparse multiply for J, M inverse, and J transpose. Note that the Rigid Body already adjusts the inverse mass and rotational inertia to be zero for pinned bodies to effectively filter (zero out) the velocity update of pinned bodies.
3. Implement the projected Guass-Seidel complementarity solver. Use the number of iterations specified in the iterations parameter. As described in class [Erleben 2007], you will want to clamp your Lagrange multiplier lambda_i to the current bounds inside your inner loop, and likewise update bounds as necessary (i.e., update bounds on the friction force tangential multiplier based on an updated normal force Lagrange multiplier).
Note that the inner loop of your iterative solver needs to be very efficient if your implementation is to work on large systems. Some of the more delicate test systems will require a larger number of iterations, smaller step sizes (e.g., higher nubmers of substeps), or perhaps some tuning of friction. You may also notice that your iterative solver converges slowly if you index your contacts in the order in which they are found. Consider randomizing the order of the inner loop in your Gauss-Seidel solver to help with this.
Demonstration movies, Novel Scene and Interesting Movie (1 mark)

Large stacks of objects will be challenging for your solver. To demonstrate the success of your implementation and how well it scales, create videos for the following test images: delicate, tower100, wallWideSparseHigh, wallWideDenseHigh. Enable drawing of center of mass positions to show how well each system comes to rest (the circles will turn solid blue when the kinetic energy falls below a threshold). Also use the mouse to perturb the system once stable to demonstrate that it isn't glue together (i.e., no bilateral constraints). Likewise, use a block transparency setting of 0.5 to help show that your Bodies are not interpenetrating.
Record a sequence of something interesting or amusing. This could be something within the existing functionality of the specification and provided code, or some additional feature you add to your code. Be creative! If you make something beautiful, consider rendering at 720p and deactivating the text overlay.

Note that you will want to turn on quality rendering features of your graphics card to reduce aliasing in your saved images. This happens when you have large images loaded. A more graceful solution would be to generate textures and mip maps for each rigid body, however, the provided code simply draws the blocks with display lists.

Challenges (optional)

Tetris factory. How many Tetris pieces can you drop before the assignment deadline? Can you fill the container? Note that complete rows (if they ever happen) do not vanish as in real game! Post snapshots to webCT to share your results with other students.

Shake a bag of nuts. Alternatively, consider using the nuts test image and add code to shake the container. Can you make the largest nuts come to the top just like in a real bag of mixed nuts? How long does it take? If the container is body zero, you can make it move by setting its velocity directly, e.g.,
```
	bodies.get(0).omega = alpha*Math.cos(simulationTime*beta);
```
How well does your simulation scale? Can you simulate thousands of bodies and thousands of contacts? Add code to report on collision solve times based on the number of contacts:
```
        System.out.println(contacts.size() + ", " + collisionSolveTime + ";" );
```
Cut and paste into a matrix A in matlab, and plot your results. Does your collision processing code scale linearly?
```
        plot(A(:,1),A(:,2),'.', 'MarkerSize', 1);
        xlabel('number of contacts');
        ylabel('PGS solve time');
```
Post your graph as a png to webCT to share, or use your graph as a test image in your simulation!
Other things to try (optional)

Broad phase collision detection. When you have large numbers of small bodies, reduce the n squared body body tests, for instance using spatial hashing.

Adaptive time stepping. Use the velocity of the fastest moving block to determine a maximum step size that prevents blocks from flying through one another. Alternatively, use continuous collision detection to prevent missing collisions.

Shock Propagation to improve stacking. See [Erleben 2007] and [Gundelman 2003] for more information on this an other optimizations.

Constraint stabilization. The LCP velocity update can also be used to update positions so that bodies are not interpenetrating. If you have errors in your velocity solve, your objects may eventually seep into one another! See [Cline Pai 2003]

Warm start for faster convergence with your projected Gauss Seidel solver. You can exploit temporal coherence of the contact forces by reusing values from your previous time step, but this involves figuring out a correspondence between contacts at different time steps.
Other forces, springs, damping, or likewise, constraint forces for joints.

References

K. Erleben, Velocity-based shock propagation for multibody dynamics animation, ACM Trans. Graph. 26, 2, Jun. 2007.

M. B. Cline and D. K. Pai, ``Post-Stabilization for Rigid Body Simulation with Contact and Constraints,'' in proceedings of the IEEE Intl. Conf. on Robotics and Autom., 2003.

Witkin, A., and Baraff, D., Eds. 2001. Physically Based Modeling: Principles and Practice. Course Notes. ACM SIGGRAPH '01. Rigid Body Simulation (slides) Constrained Dynamics (slides)

Eran Guendelman, Robert Bridson, Ronald P. Fedkiw, Nonconvex Rigid Bodies With Stacking, ACM Transactions on Graphics, 22(3), July 2003, pp. 871-878.

Cornell CS2643 rigid body contact assignment by Doug James.

Finished?