COMP259: Homework 3

Simple Constrained Dynamics

Problem A: A Bead on a Wire

Simulate the motion of a particle and visualize the result. The particle is defined in a Cartesian coordinate system p=(x,y), and it is constrained by a unit circle centered at the origin. Use the initial position (0, 1), the initial velocity (0, 0), and gravitational force (0,-mg). Constants are give as m=1Kg, and g=9.8m/sec2.

Theory

To enforce the constraint holding the bead on the wire I used the Lagrangian Dynamics formulation presented in the SIGGRAPH'97 course notes on Physically Based Modeling by David Baraf and Andrew Witkin. Lagrangian Dynamics involves representing the object in a coordinate system that must necessarily satisfy the constraints. In the case of the bead on a wire, we represent the bead's state as simply an angle theta and an angular velocity about the circle. The difficulty in the problem then becomes how to convert forces that are represented in standard Euclidean coordinates into accelerations in the angle theta. The course notes provide the derivation of the equation used to do this. To apply this general framework to our specific example we need to compute the Jacobian of the coordinate transformation and its time derivative. But once this has been worked out on paper, I found that many of the terms of the equation cancel, leaving the simple formula:

angular_acceleration  =( f_y*cos(theta) - f_x*sin(theta)) / (mass * radius)

where f_x and f_y are the components of the external force in 2D Euclidean coordinates, theta is the current angle of the particle, mass is the particles mass, and radius is the radius of the wire circle. This angular acceleration can then be integrated using any integration method to produce the correct motion of the bead on the wire, in response to the external force.
 
 

Implementation

Particle Dynamics

I implemented the Lagrangian Dynamics, described above, in a class called PBMCircleParticle. It uses the integrator framework developed for Homework 1 to integrate the dynamics for the angle theta. The function Update(timeStep) packs the state data into a configuration vector for the ODE solver and invokes the solver to obtain the state for the next time step. The function ODEDerivative() is used as a callback by the ODE solver to obtain the derivative of the system. This function is where the equation described above, for computing the angular acceleration from a Euclidean force is implemented.

Source: PBMCircleParticle.h, PBMCircleParticle.cpp
 

Main Program

The main renders the scene, calls the routine to update the particle and allows user interaction. The user can change the radius of the circle (PgUP & PgDwn), change the integration method (key "i"), and add a directional force  to the system in the direction of the mouse (key "f"). The force added to the system can be scales using the Up and Down arrow keys.

Source: main.cpp