GSoC, Week 4

Things went pretty well this week.  Last weekend, I implemented a Dyad class.  We decided to do this for a number of reasons; but largely it was to allow the end-user the freedom to express a body’s rotational inertia as sums of components in different frames.  More information about dyads can be found here:, or in Kane’s 1985 book: Dynamics, Theory and Applications.  A dyad (or dyadic, I guess?) represents the juxtaposition of two vectors.  My experience with dyadics has involved treating them as a matrix, basically; but implementing this class showed me some cool things they can do. I also made a pretty simple convenience function for creating an inertia dyadic by supplying 3 or 6 scalar values and a frame; this will allow users to treat the inertia dyadic like they would treat the inertia tensor; but they could still write the more complicated Dyad expressions out if they wanted.  Also, the outer product was implemented for Vector, allowing for easy creation of Dyads.

So, Dyad was the last main class which is used to describe physical parameters/quantities/etc. After that, RigidBody and Particle classes were made.  Luke and I discussed how we implement Kane’s Method when doing it by hand, and realized that basically we do all of the kinematics, then just write down tables for storing the partial velocities and forces.  We both came to the realization that the terms “particle” and “rigid body” were nothing more than associating a mass/inertia with a point/frame, and that forces are really just a Vector associated with a Point or ReferenceFrame.  So, we decided to make RigidBody and Particle container classes effectively, where they both have attributes (and getters/setters for sanitizing input) for all their relevant information.  Then, I decided that storage of forces didn’t have to be anything more than a list of tuples, in the form (point/frame, force/torque).  I think we’ll make some more convenience functions here, probably relating to gravitational forces.

Finally, Kane’s Method.  Luke has done a lot of thinking about this, so a lot of credit goes to him for the successful implementation of this.  When one does Kane’s Method by hand (or at least when I do it), I follow a series of steps; something along the lines of:

  1. Set up generalized coordinates and speeds
  2. Do kinematics
  3. List forces
  4. Form relevant partial velocities
  5. Calculate Fr, Fr*
  6. Rewrite in desired form
So, I’ve basically written a class which stores the relevant information as you go through these steps. Here’s a little example:
KM = Kane(N)
KM.gen_speeds([u1, u2, u3, u4, u5, u6])
KM.dependent_speeds([u1, u4, u6], conl)
fr = KM.form_fr(FL)
frstar = KM.form_frstar(BL)

First, you create the object specifying the inertia frame.  Then you supply the generalized speeds and kinematic differential equations, which relate all of the q’s (gen. coords.) to u’s (gen. speeds) by means of qdot = f(u).  This is stored in a dictionary, and I need to figure out a better/easier way for the user to generate these for bodies oriented by Euler Angles (perhaps another convenience function?).  Fr is the generalized active force, and is formed from the list of forces/points and their partial velocities, Fr* is the generalized inertia force, and is formed from the list of particles/bodies and their partial velocities.

Something that Luke has thought about for a while is dealing with nonholonomic systems; these are systems with velocities constraints (example: ideal skate, where it can move forward and backward and it can yaw, but never move left and right).  Kane describes it in his 1985 book, but it’s a little unclear; the way Luke has described is the same thing, but a little easier conceptually.  This will get described in a lot more detail in the documentation for PyDy.  One other cool thing which I wanted to implement (so I did) was dealing with systems with time varying mass.  Newton’s second equation (f = m a) is really f = d/dt (m v); or force = time derivative of linear momentum.  This also extends to angular momentum/torques.  I make sure to check the time derivative of the mass or inertia, and if it is zero, do things the “normal” way; if it is time varying, I calculate the momenta and take the time derivatives thereof.  So, my code should work with systems with non-constant mass/inertia.  Finally, one more cool thing.  Sometimes, a system will be defined by parameters which are time-varying, but not dynamics (ie user-specified speeds or positions or such).  These will be defined as DynamicSymbols instead of Symbols.  Then, when there are DynamicSymbols which are not generalized speeds or coordinates, the code will find and identify them.  The plan is for the SciPy output code to have empty functions for these parameters for the user to then fill in.  One thing I haven’t thought about (mainly because I have never done any examples) is flexible dynamic systems; I’d really like for my code to work with flexible systems, but I’ll need to try a few first.  Hopefully it will transfer nicely into the framework I have built so nothing new will have to be implemented.

So basically, Kane’s Method has been written.  The functionality that was supposed to be complete this week was supposed to include SciPy code generation, but I’m a little behind due to my late start.  I feel like I’ve made a good amount of progress though…… The schedule for next week has creation of equations of motion unittests, so I’ll probably work out a bunch of examples (or take them from Kane’s book) and implement them as tests, in addition to implementing SciPy code generation.  The week after is integration into SymPy, so I’ll probably start talking to Aaron this week about that.  Then it’s examples, documentation and LaTeX output.  Next week I’ll have some examples of generated EoM to show.