GSoC, Week 2

My blog post is a little late this week…

Good progress was made this week, with lots of code added to the kinematics file, starting a tests file, and starting the inertias file.  I think most of the kinematics classes are done…. but more tests need to be implemented, and then the code corrected if there are any errors.  It’s looking like the inertias file will be pretty easy to fill out, as those classes are extending a lot of what has already been written.  The next code that will be a challenge to figure out is the algorithm implementations (e.g. Kane’s Method, Lagrange’s Method, etc.).  In forming the equations of motion, expressions can get very big very fast, so we will have to take care here.

I think the best piece of news was from yesterday.  Throughout the week, I had been trying to figure out how to be able to take the partial derivative of a SymPy expression with respect to a generalized speed, without using any substituions (that the user would see, at least).  The problem was that SymPy’s diff can only take in symbols, and for an undefined function (say, x(t)), you get a Derivative object returned (say, D(x(t), t)).
At first I tried extending Symbol, unsuccessfully, to return a new symbol when its derivative was taken.  This didn’t work once my extended Symbol was inside a Add or Mul, as the new methods I had written weren’t called.  I posted this to the mailing list (http://groups.google.com/group/sympy/browse_thread/thread/b9c8e3fc3d3d379b), and after some discussion, someone suggested that I add Symbol(‘t’) to the free_symbols property of my extended Symbol class.  I also got the name for the class, DynamicSymbol, from this discussion.  Basically, now when you take the time derivative of a DynamicSymbol object, even if it is within other SymPy objects, it’s _eval_derivative() method is called, and now it will return what you want.  Here’s a little example:

In [1]: run kinematics.py

In [2]: x = TVS(‘x’)

In [3]: y = TVS(‘y’)

In [4]: t = Symbol(‘t’)

In [5]: diff(2 * x**2 + 4 + y,t)
Out[5]: 4*x*xdt + ydt

In [7]: diff(diff(2 * x**2 + 4 + y,t),TVS(‘xdt’))
Out[7]: 4*x

This is pretty exciting, as now things will work fairly intuitively, when it comes to derivatives and partial derivatives for generalized coordinates and speeds.  Also, if you have time-varying symbols (maybe some specified position), it will now be identified as such, and when the final equations of motion are formed with the output code, it can be treated as a variable defined by a user-specified function.

Next is filling out the inertia classes, adding more tests and ensuring correct functionality, and starting on the algorithm implementations.

GSoC, Week 1

Well, this week was the start of Google Summer of Code 2011.  Unfortunately, the quarter is still not over for me, so I didn’t make as much progress as I would like.  But, it looks like some of the questions about rewriting the Vector class have been answered.

When doing dynamics problems by hand, basis vectors are typically written as:

  or  

Where e is the reference frame within which the basis vectors are “standard”.  One thing that is immediately obvious is that the numbering starts from 1, unlike most programming languages, which start indices at 0.  PyDy is going to be printing out equations to LaTeX, for use in publication; this would most certainly use 1 indexing.  This leads to some conflicting things, where e[0] in the code is in LaTeX.  This is something I definitely did not want to happen.  So e.x, e.y, and e.z will be how the basis vectors are accessed, removing the indexing question.  There will be a LaTeX printing option to determine 123, xyz, or ijk indices.  But, there will never be 0 as an index (that the end-user will see).

My decision to keep the basis vectors in the ReferenceFrame class, versus generating them upon some function call, ultimately came down to the same issue.  If you named a set of basis vectors b1,b2,b3 , but they were in frame “Foo”, printing them would give foo1>, foo2>, foo3>.  This again, would be unacceptable.  This situation is also fairly likely to come up, as one usually wants to take shortcuts (I know I have when using Autolev).  It also ties into what I guess is my last point for this topic: readability for non-python programmers.

The professor who will teach the graduate multibody dynamics course next year said that he would consider teaching the students PyDy next winter quarter.  This is going to mean that the code should be readable by people who aren’t familiar with Python.  Here, by code, I mean the scripts that people will write to generate their equations of motion.  Anyways, readability for non-Python users, while still trying to be Pythonic, will be an important consideration this summer.

Next is finishing up the Vector and ReferenceFrame classes.  I’m trying to reuse as much code as possible, and definitely not lose any functionality.  It seems to be going OK so far and hopefully it will stay that way.

PyDy and other things

This weekend is Maker Faire, so the lab is going there.  One thing I’ll have is some info on some MEMS sensors from a class project last quarter.  We collected info on the Analog Devices ADXL 345 and Invensense ITG-3200 (tried to use a Honeywell HMC5843, but we couldn’t read it over I2C).  Anyways, I got permission from my group members to post our report online, so here it is.

Report: MAE276 Final Paper

For PyDy/SymPy & GSoC there is also some news.  It looks like the previous separation between UnitVector and Vector classes is going away; in addition they will no longer extend SymPy’s Basic or Expr classes.  There will now be only one Vector class.

It will store a list of lists; the inner list will have the 3 vector measure numbers for a frame (something like [([1],[2],[4]),’b’] ), where the measure numbers will use the SymPy Matrix class (and can have symbols in them).  The outer list will hold each of these lists for each frame.  So the data stored will be something like [[([1],[2],[4]),’b’],[([3*x],[1],[0]),’c’]]… or something like that.  The current plan is to write our own operators for addition, subtraction, scalar multiplication & division, and dot and cross products.  By not using the SymPy classes, we can better control how the Vector class will behave.  One important behavior is to have every operator (except for dot product) take in and return a vector, and not have any confusion between SymPy objects and Vectors when using operators.

There is also some debate as to whether the ReferenceFrame class should store its own basis vectors, or they should be returned upon initialization (see the sympy list post: Question regarding vectors.  I understand the argument to keep them with their frame and the awkwardness that will exist when creating the basis vectors upon ReferenceFrame initialization.  But, I feel that the Vector class already stores information about each frame (in its inner lists), and question whether basis vectors are a property of ReferenceFrames, or a frame is a property of a set of basis vectors.

I guess the orthonormal properties of a set of 3 basis vectors exist independently of a particular reference frame.  I guess what defines a frame?  I think of it as a rotation (and related time derivatives) from on set of basis vectors to another.  So it is defined by both its basis vectors and a rotation.  I still feel like its is preferable to keep the basis vectors outside of the ReferenceFrame class, but not sure my position has a strong enough argument (but I also don’t think it is too weak of an argument).  Hopefully starting to discuss more of the ReferenceFrame class methods will clear this up.

Equations of Motion Example

Here, I have derived by hand the equations of motion for a rolling disc.  The scanned document at the bottom shows the steps I went through.

This is to show what PyDy (my Google Summer of Code 2011 project, under SymPy) will be doing; by using SymPy’s symbolic manipulation abilities, one can use the functions in PyDy to derive equations of motion symbolically.  A simple system like the rolling disc can be done by hand easily but for larger, multibody systems you will need some computational assistance.

PyDy is primarily built around Kane’s Method ¹ , but will be able to work with Lagrange’s Method or Newton-Euler.

Edit: PDF removed for now until I fix an error….

1. Kane, T. R. and Levinson, D. A. (1985)  Dynamics: Theory and Applications. MacGraw-Hill Series in Mechanical Engineering. MacGraw-Hill Book Company, NewYork.