348C Assignment #2: Constrained Dynamics (inextensibility
Professor: Doug James
Due date: Thursday, Oct 26, 2017 (midnight)
Inextensible yarn curves from [Kaldor et al.
In this second assignment, you will extend your assignment #1 code
to support hard holonomic (position-based) constraints.
Specifically, you will use constraints to implement deformable 2D
strands with inextensible (fixed-length) line segments. Why? Well,
stretchy spaghetti looks weird, and (as you will see later)
exploiting these constraints can speed up other calculations, e.g.,
of bending forces. You will use your previous collision detection
and resolution implementation, and should be able to re-run your
prior examples with better length preservation.
Background: As discussed in class, inextensibility
constraints will be added for each edge to maintain a fixed length
equal to its rest length. For example, we will use constraints to
enforce that the squared-length of an edge's inter-particle distance
is constant. The formulae for constraints, Jacobian, and constraint
matrix system were covered in class.
Implicit constraint direction: There are many ways to
discretize and enforce constraints. In this assignment, you will use
an implicit velocity-level formulation that seeks to apply an
impulse to our predictor velocities (from symplectic Euler) such
that the constraint will be enforced at the end of the
timestep. The details of this impulse computation have been
described in class (Tues Oct 17), and are elaborated on further
below; for more details on ICD for cloth simulation, see [Goldenthal et
Combining collisions with hard constraints: As they say, when
an immovable object hits an immovable wall, something has to give.
Unless you simultaneously resolve both the hard constraints and the
collision constraints, they will not both be satisfied.
Nevertheless, we can alternate between the two classes of
constraints, to improve their satisfaction. In your
implementation, you will (a) solve all inextensibility constraints
simultaneously using a matrix solve (details below), and (b) perform
the velocity update due to the inextensibility constraint just
before the main collision processing loop, but after the velocities
have been updated with other forces (gravity, springs, penalty
forces). The inextensibility constraints will be satisfied for
all edges, unless they are adjacent to vertices involved in contact
impulses. Once you get this working, you can further iterate on both
inextensibility and collision constraints (using velocity-level
impulse corrections for both), however, this may not be necessary
due to our ICD formulation.
How to demonstrate the merits of your simulator:
Matrix Solver Details: To greatly simplify your your
implementation, you need not support (1) edges with only one pinned
vertex, or (2) branching edge structures. In other words, you only
scan to find strands in free-flight motion. For each strand, the
constraint matrix is a positive definite, tridiagonal matrix, which
can be solved efficiently using a linear-time solver such as a
classical tridiagonal matrix algorithm (TDMA), such as the Thomas
Your implementation should identify the vertices in each separate
strand (ideally once, at the beginning of the simulation), then
assemble and solve each linear system using a TDMA, the resulting
impulses are then used to update the particle velocities prior to
impulse-based collision processing. (Due to repeated reuse, you
should be able to optimize your code to avoid creating/destroying
arrays, and simply track strands. Note that the associated
constraint matrix and RHS are position dependent, and therefore
likely to change at each timestep.)
- Simulate a single very long strand of several dozen particles
to the one on the A1 webpage).
- Report the maximum relative length error (|delta
Length|/|length|) for all constrained edges for all time, e.g.,
as the strand falls on the bottom of the container.
- After you've done (1) & (2), re-run the spaghetti factory
to demonstrate the effect of inextensibility. There is no need
to run for a very long time, however you must identify all
strands and enforce constraints on each one. Since meatballs
create a constraint loop, they are not handled by the TDMA, and
therefore you may simply disable the MeatballFactory for this
ADDITIONAL THINGS TO TRY (Optional): Only if you have time,
or for a future (final) project, here are some logical next things
Hand-in using Canvas:
- Constraints: Support pinned-vertex constraints in your
constraint solver. You can do this by adding additional
constraints, in more than one way.
- Constraints: Support branching or looping structures,
e.g., meatballs. This will require a better matrix solver, which
- Numerical Linear Algebra: Try a different matrix
solver, such as Cholesky factorization (direct solver) or the
preconditioned conjugate gradient (PCG) method (iterative
- Numerical Linear Algebra: Instead of a direct matrix
solver, try an iterative relaxation scheme (such as Jacobi or
Gauss-Seidel) that essentially relaxes just one constraint at a
time. Iterative until you reach a fixed number of iterations, or
achieve a desired relative-error tolerance.
- Constraint Drift: How well are your inextensibility
constraints being satisfied? Do you notice a difference if you
enable/disable the edge spring forces from A1?
- Modeling: Try isometric bending forces that exploit the
fact that edge lengths are constant when deriving the bending
forces. Note this 2D optimization is analogous to the 3D cloth
bending forces for which isometric deformations (with
inextensible triangle edges) lead to fast quadratic bending
et al. 2006]. Note that when the inextensibility
constraint is satisfied, there should be no difference in your
bending forces except that they are much faster. One further
benefit is that an implicit integrator is more easily
implemented given the simpler derivative formula.
- Modeling: Simulate a furry creature. You can ignore
collisions for real-time speed, but you will need to fix the end
position and orientation of each hair.
Start early. Ask questions. Have
- Document: A very brief
written document (in txt or PDF) that details what you did, your
findings, any problems, and who you discussed the assignment
- Include any comparisons or analysis/results that demonstrate
the functioning behavior of your implementation---primarily
the "largest relative length error" information.
- Include the names of people that you discussed the
assignment details with.
- Include citations of external resources (books, webpages,
etc.) that you used in your R&D.
- Code: Documented
software implementation derived from the Java starter code, or
based on an analogous framework in another language, e.g., C++.
- Results: Your
creative simulation artifacts, videos, images,
- We may run your software to evaluate your implementation,
and/or ask you to provide a demo, wowever please
submit video that demonstrates the
- Submit as a group: If
you are working with a partner, be sure to form a two-person
group, and submit your zip file as a group.
- Late submission: Keep in mind this
is a short assignment which should take less time than other
assignments. Please submit on time so that we can continue to
the next assignment.
On collaboration and academic
are allowed to collaborate on the assignments to the extent of
formulating ideas as a group, and derivation of physical equations.
However, you must conduct your programming and write up completely
on your own (or with your partner), and understand what you are
writing. Please also list the
names of everyone that you discussed the assignment details
are expected to maintain the utmost level of academic integrity in
Andrew Witkin, Physically Based Modeling, Constrained
Dynamics, Online SIGGRAPH 2001 Course Notes, 2001.
- Rony Goldenthal,
David Harmon, Raanan Fattal, Michel Bercovier, Eitan Grinspun, Efficient Simulation of
Inextensible Cloth, ACM Transactions on Graphics,
26(3), July 2007, pp. 49:1-49:7. [ACM Digital Library link]
M. Kaldor, Doug L. James, Steve Marschner, Simulating Knitted Cloth at the
Yarn Level, ACM Transactions on Graphics,
27(3), August 2008, pp. 65:1-65:9. (uses inextensibility
constraints on yarn spline curves)
- Miklós Bergou, Max Wardetzky, David Harmon, Denis Zorin, Eitan
Quadratic Bending Model for Inextensible Surfaces,
Fourth Eurographics Symposium on Geometry Processing, pp.
227-230, June 2006. (Project
- U.M. Ascher and
L.R. Petzold, Computer Methods for Ordinary
Differential Equations and Differential-Algebraic Equations,
SIAM. (good introductory book on ODEs & DAEs)
- E. Hairer and G.
Wanner, Solving Ordinary Differential
Equations II: Stiff and Differential-Algebraic Problems,
2nd edition, Springer, 1996. (See
Chapter VII.(1-2) Differential-Algebraic Equations of Higher
Index). (excellent but advanced text on DAEs)