Professor: Doug James

Due date: Friday, March 2 (before midnight).

Videos: Please submit videos (described below) that demonstrate the requested (and any additional) functionality of your system. Details on video generation are here.

In this assignment you will add implicit integration to your mass-spring particle system, modifying your first assignment code appropriately. For simplicity, you are only required to integrate a single undamped hair model implicitly. Using your implementation, you will be able to simulate a long strand of stiff hair falling on the ground and colliding with the walls of the cell, as well as interact with it interactively (include a final video of this).

Assignments Steps: Here are the steps you need to address:

- Procedural hair generation:
Automate the construction of a hair primitive parameterized by the
number of particles, N--so you don't have to draw it by hand each
time. Make it so that your simulator can easily create a long
hair with a moderate N value, e.g., N=50-500. With stiff spring
forces, you will notice that your explicit integrator becomes quite
inefficient and/or unstable for such models.

- Backward Euler Integrator:
You will integrate the hair's stiff extension and bending spring forces
using backward Euler as in the Baraff course notes (Implicit
Methods for Differential Equations). The equation that you
will need to set-up and solve, hereafter
called A*dv=b, is
given in (4-6)--note that f in
(4-6) is actually the particle acceleration vector. For simplicity we
will ignore spring damping forces, and
therefore the acceleration function f
will only depend on the particle positions, x. Therefore you will only need to
compute the function, f0, and
it's sparse Jacobian matrix, df/dx.
How and where you accumulate the Jacobian is something you should think
about.

- Spring
force derivatives: The first step is to add sparse derivative support to
the Force objects used to model a single hair (
**SpringForce2Particle, SpringForceBending)**. You can modify the Particle object to accumulate sparse force Jacobians as well forces (.f), keeping in mind that the Jacobian is a sparse matrix with nonzero 2x2 blocks. You could store/access this sparse block row vector in a variety of ways, including just using a HashMap<Integer,Matrix2d> accumulator--keep in mind that the Jacobians are symmetric. (If you go this route, you'll probably want to make a helper object that handles get/set/acc/multiply operations involving the indexed Matrix2d information and global double[] vectors.) You could then modify the Force interface to optionally support apply_forceAndJacobian() to accumulate both force and Jacobian information. Alternately, you can accumulate the sparse derivatives in a generic sparse matrix container from an external library, e.g., MTJ. It's really up to you. - Conjugate
Gradients: Once you can compute the quantities in equation
(4-6), then you can solve the system, A*dv=b. Since A is
symmetric and positive definite we will use the Conjugate Gradient (CG) iterative method.

- Note that it is not required to formally assemble the matrices df/dx, or A=(I - h^2*df/dx) in a sparse matrix container in order to solve the matrix equation: you only need to evaluate b=h(f0 + h*df/dx*v0), then provide CG with a way to evaluate matrix-vector products Az, for some vector z. For example, here is a simple java interface for a matrix-vector multiply function that could be implemented by some object (that uses your ParticleSystem to appropriately accumulate force/Jacobian information, and then repeatedly/quickly multiply by it):
- public
interface MatrixMultiplyInPlace {

/**

* Multiplies a double[] and puts the result in a preallocated double[].

* @param x Vector to be multiplied by matrix.

* @param Ax Matrix-vector multiply result.

*/

public void multiply(double[] x, double[] Ax);

} - Note that you can implement
this interface to multiply by the total spring force df/dx, as well as to implement A by composition.

- Feel free to use any implementation of the CG algorithm that you would like, however it is quite easy and instructive to "roll your own": simply type in the implementation from Shewchuk's paper (page 50... it's not very painful ;). To manipulate double[] vectors you'll need some kind of simple java BLAS (add/acc/dot/axpby/twoNorm), and here is a class (djames.matrix.DoubleVector) which has functions you can borrow.
- Accuracy
is something that you can adjust to make your simulation
go faster, but be careful since for large relative errors you'll start
to see artifacts. For interactive applications, people often only take
a fixed number of iterations per time-step to avoid big delays, but
this isn't guaranteed to be accurate. Another problem is that the
A matrix may not be
positive for large time-steps (doh!), and so CG may fail to converge.

- Preconditioned CG (PCG): In practice, CG can converge a lot faster if you use almost any reasonable preconditioner (see Shewchuk, page 51--it's only ~two extra lines which isn't too painful). For example, just using the inverted block diagonal can be a big help. Try to use a simple preconditioner, but it is not required. You can have it implement MatrixMultiplyInPlace for use in PCG (reverting to CG if the preconditioner is null).
- Additional
external forces: You can add additional forces to the system,
but you need not integrate them implicitly (phew!). For example, user
interactions can be modeled using an explicitly integrated spring
force, as can gravitational forces, e.g., using forward Euler. Just remember which
Force objects are explicit and implicit in your "advanceTime()" method.

- Collisions
and other constraints:
Collisions and other constraints, e.g., Pin constraints, will be
handled using a velocity-level filter, as (will be) described in class.
You should be able to pin your models, as well as have them collide
with the walls of the unit computational cell.

- Create Artifacts: See if you can do something interesting. For example, can you model a head of hair? a furry beast? some long grass or reeds? Let your imagination run wild.
- Videos:
Please
include videos illustrating the behavior of your system for a single
long strand of hair, and if
possible its interactive performance for a moderately sized system,
e.g., N=100----how big a strand can you do at 10 FPS?

Hand-in using CMS: Please submit a short write-up (detailing what you did, your findings, and who you discussed the assignment with, etc.), as well as your Java implementation, and any artifacts, videos, etc. Submit videos in a compressed format.

Enjoy!!

Copyright Doug James, February 2007.