Permalink
Find file Copy path
Fetching contributors…
Cannot retrieve contributors at this time
699 lines (662 sloc) 45.7 KB

PyDy for Autolev Users

Introduction

Autolev (now defunct) is a domain specific programming language which is used for symbolic multibody dynamics. The SymPy mechanics module now has enough power and functionality to be a fully featured symbolic dynamics module. The PyDy package extends the SymPy output to the numerical domain for simulation, analyses and visualization. Autolev and PyDy have a lot in common but there are also many differences between them. This page shall expand upon their differences. It is meant to be a go-to reference for Autolev users who want to transition to PyDy.

It would be helpful to have a basic idea about scientific computing with Python before going over this page. If you aren’t familiar with scientific computing with Python there are many sources to learn. You can start with the Python programming language itself, with the canonical source being the Python Documentation. The SciPy Lectures are a great intro to scientific computing with Python. Finally, to learn about how to do symbolic manipulation with SymPy, check out the SymPy Documentation, especially the tutorial.

Some Key Differences

Autolev PyDy
Autolev is a domain specific programming language designed to perform multibody dynamics. Since it is a language of its own, it has a very rigid language specification. It predefines and computes many things based on the input code. Its statements are a lot cleaner as a result of this.
PyDy is a library written in the general purpose language Python. Although Autolev's code is more compact, PyDy(by virtue of being an add on to Python) is more flexible. The users have more control over what they can do. For example, one can create a class in their code for let's say a type of rigibodies with with common properties.
Autolev generates Matlab, C, or Fortan code from a small set of symbolic mathematics.
PyDy generates numerical Python, C or Octave/Matlab code from a large set of symbolic mathematics created with SymPy. It also builds on the popular scientific Python stack such as NumPy, SciPy, IPython, matplotlib, Cython and Theano.
Autolev uses 1 (one) based indexing. The initial element of a sequence is found using a[1].
Python uses 0 (zero) based indexing. The initial element of a sequence is found using a[0].
Autolev is case insensitive.
PyDy code being Python code is case sensitive.
One can define their own commands in Autolev by making .R and .A files which can be used in their programs.
PyDy code is Python code, so one can define functions in their code. This is a lot more convenient.
Autolev is proprietary.
PyDy is open source.

Rough Autolev-PyDy Equivalents

The tables below give rough equivalents for some common Autolev expressions. These are not exact equivalents, but rather should be taken as hints to get you going in the right direction. For more detail read the built-in documentation on SymPy vectors , SymPy mechanics and PyDy .

In the tables below, it is assumed that you have executed the following commands in Python:

import sympy.physics.mechanics as me
import sympy as sm

Mathematical Equivalents

Autolev PyDy Notes
Constants A, B
a, b = sm.symbols (‘a b’, real=True)
Note that the names of the symbols can be different from the names of the variables they are assigned to. We can define a, b = symbols(‘b a’) but it is good practice to follow the convention.
Constants C+
c = sm.symbols(‘c’, real=True, nonnegative =True)
Refer to SymPy assumptions for more information. modules/core
Constants D-
d = sm.symbols(‘d’, real=True, nonpositive=True)

Constants K{4}
k1, k2, k3, k4 = sm.symbols('k1 k2 k3 k4', real=True)

Constants a{2:4}
a2, a3, a4 = sm.symbols('a2 a3 a4' , real=True)

Constants b{1:2, 1:2}
b11, b12, b21, b22 = sm.symbols('b11, b12, b21, b22', real=True)

Specified Phi
Phi = me.dynamicsymbols(‘Phi ')

Specified F = k1*x + k2*theta + k3*x' + k4*theta'
F = me.dynamicsymbols('F')
F = k1*x + k2*theta + k3*xd + k4*thetad

Variables q, s
q, s = me.dynamicsymbols (q, s)

Variables x’’

x = me.dynamicsymbols(‘x’)

xd = me.dynamicsymbols (‘x’, 1)

xd2 = me.dynamicsymbols(‘x’, 2)

 
Variables y{2}’

y1 = me.dynamicsymbols (‘y1’)

y2 = me.dynamicsymbols (‘y2’)

y1d = me.dynamicsymbols(‘y1’ , 1)

y2d = me.dynamicsymbols(‘y2 , 1)

 
MotionVariables u{2}

u1 = me.dynamicsymbols(‘u1’)

u2 = me.dynamicsymbols(‘u2’)

SymPy doesn’t differentiate between variables, motionvariables and specifieds during declaration. Instead, it takes different lists of these as parameters in objects like the KanesMethod.
Imaginary j

j = sm.symbols(‘j’)

j = I

I is a sympy object which stands for the imaginary unit. One can define complex numbers using it.

z = x + I*y

where x, y and z are symbols.

Tina = 2*pi
Tina = 2*sm.pi
If one wants to use numerical constants instead of symbolic ones they can import them from mpmath.
abs(x)^3 + sin(x)^2 + acos(x)
sm.abs(x)**3 + sm.sin(x)**2+ + sm.acos(x)

E = (x+2*y)^2 + 3*(7+x)*(x+y)

Expand(E, n:m)

Factor(E, x)

Coef(y, x)

Replace(y, sin(x)=3)

E = (x+2*y)**2 + 3*(7+x)*(x+y)

sm.expand(E)

sm.factor(E)

y.coeff(x)

y = y.xreplace ({sm.sin(x): 3})

For more information refer to simplification.

Dy = D(E, y)

Dt = Dt(E)

sm.diff(E, y)

sm.diff(E, t) where t = me.dynamicsymbols._t

Works if the expression is made up of dynamicsymbols.

For more information refer to calculus.
TY = Taylor(x*cos(x), 0:7, x = 0) ty = taylor(x*sm.cos(x) ,0 , 7)

Execute

from sympy.mpmath.import *

For more information refer to mpmath/calculus.

F = Evaluate(E, x=a, y=2)

E.subs([(x, a), (y, 2)])

To get floating point numbers from numerical expressions use evalf

E.evalf((a+sm.pi).subs ({a:3}))

 
P = Polynomial([a, b, c], x)
p = sm.Poly(a*x**2 + b*x + c)
For more information refer to modules/polys.
Roots(Polynomial([a, b, c], x), x, 2) sm.solve( sm.Poly(a*x**2 + b*x + c))

For more information refer to solvers.

For numerical computation related to polynomials and roots refer to mpmath/calculus.

Solve(A, x1, x2)

where A is an augmented matrix that represents the linear equations and x1, x2 are the variables to solve for.

sm.linsolve(A, (x1, x2))

where A is an augmented matrix

For more information refer to solvers/solveset.

RowMatrix = [1, 2, 3, 4]

ColMatrix = [1; 2; 3; 4]

MO = [a, b; c, 0]

MO[2, 2] := d

A + B*C

Cols(A)

Rows(A)

Det(A)

Element(A, 2, 3)

Inv(A)

Trace(A)

Transpose(A)

Eig(A)

row_matrix = Matrix([[1],[2], [3],[4]])

col_matrix = Matrix([1, 2, 3, 4])

MO = Matrix([[a, b], [c, 0]])

MO[1, 1] = d

A + B*C

A.shape(0)

A.shape(1)

M.det()

M[2, 3]

M**-1

trace(A)

A.T

A.eigenvals()

For more information refer to matrices.

Physical Equivalents

Autolev PyDy Notes

Bodies A

Declares A, its masscenter Ao, and orthonormal vectors A1>, A2> and A3> fixed in A.

m = symbol(‘m’)

Ao = symbols(‘Ao’)

Af = ReferenceFrame(‘Af’)

I = outer(Af.x, Af.x)

P = Point(‘P’)

A = RigidBody(‘A’, Ao, Af, m, (I, P))

Af.x, Af.y and Af.z are equivalent to A1>, A2> and A3>.

The 4th and 5th arguments are for the mass and inertia. These are specified after the declaration in Autolev.

One can pass in None for the parameters and use setters A.mass = _ and A.inertia = _ to set them later.

For more information refer to mechanics/masses.

Frames B
B = ReferenceFrame(‘B’)
For more information refer to physics/vector.
Newtonian N
N = ReferenceFrame(‘N’)
SymPy doesn’t specify that a frame is inertial during declaration. Many functions such as set_ang_vel() take the inertial reference frame as a parameter.
Particles C

m = symbol(‘m’)

po = Point(‘po’)

C = Particle(‘C’, po, m)

The 2nd and 3rd arguments are for the point and mass. In Autolev, these are specified after the declaration..

One can pass in None and use setters (A.point = _ and A.mass = _) to set them later.

Points P, Q

P = Point(‘P’)

Q = Point(‘Q’)

 
Mass B=mB

mB = symbols(‘mB’)

B.mass = mB

 
Inertia B,I1,I2,I3,I12 I23,I31

I = inertia(Bf, i1, i2, i3, i12, i23, i31)

B.inertia = (I, P) where B is a rigidbody, Bf is the related frame and P is the center of mass of B.

Inertia dyadics can also be formed using vector outer products.

I = outer(N.x, N.x)

For more information refer to the mechanics api.

vec> = P_O_Q>/L

vec> = u1*N1> + u2*N2>

Cross(a>, b>)

Dot(a>, b>)

Mag(v>)

Unitvec(v>)

vec = (Qo.pos_from (O))/L vec = u1*N.x + u2*N.y

cross(a, b) where a and b are vectors

dot(a, b)

v.magnitude()

v.normalize()

For more information refer to physics/vector.

P_O_Q> = LA*A1>

where O is a point

P_P_Q> = LA*A1>

where P is a particle

Q.point = O.locatenew(‘Qo’, LA*A.x)

where A is a reference frame.

Q.point = P.point.locatenew(‘Qo ’, LA*A.x)

For more information refer to the kinematics api.
V_O_N> = u3*N.1> + u4*N.2>
O.set_vel(N, u1*N.x + u2*N.y)

A_O_N> = 0>

Acceleration of point O in reference frame N.

O.set_acc(N, 0)  

W_B_N> = qB’*B3>

Angular velocity of body B in reference frame F.

B.set_ang_vel(N, qBd*Bf.z)

where Bf is the frame associated with the body B.

 

ALF_B_N> = Dt(W_B_N>, N)

Angular acceleration of body B in reference frame N.

B.set_ang_acc(N, diff(B.ang_vel_in(N))  

Force_O> = F1*N1> + F2*N2>

Torque_A> = -c*qA’*A3>

In SymPy one should have a list which contains all the forces and torques.

fL.append((O, f1*N.x + f2*N.y))

where fL is the force list.

fl.append((A, -c*qAd*A.z))

 

A_B

or

Dircos(A,B)

where A and B are reference frames

A.dcm(B)  
CM(B)
B.masscenter

Mass(A,B,C)
A.mass + B.mass + C.mass

V1pt(A,B,Bq,Q)
Q.v1pt_theory(Bq, A, B)

V2pts(A,B,P,Q)
Q.v2pt_theory(P, A, B)

A1pt(A,B,Bq,Q)
Q.a1pt_theory(Bq, A, B)

A2pts(A,B,P,Q)
Q.a2pt_theory(P, A, B)

Angvel(A,B)
B.ang_vel_in(A)

Simprot(A, B, 1, x)
B.orient(A, ‘Axis’, x, A.x)

Gravity(G*N1>)
fL.extend(gravity(g*N .x, P1, P2, ...))

Force(P/Q, v>)
fL.append((P, -1*v), (Q, v))

Torque(A/B, v>)
fL.append((A, -1*v), (B, v))

Fr()

FrStar()

(fr, frstar) = KM.kanes_equations(fL , bL)

where KM is the KanesMethod object.

For more details refer to mechanics/kane and the api.
Kindiffs(A, B ...)
KM.kindiffdict()

Momentum(option)

linear_momentum(N, B1, B2 ...)

reference frame followed by one or more bodies

angular_momentum(O, N, B1, B2 ...)

point, reference frame followed by one or more bodies

 
KE()

kinetic_energy(N, B1, B2 ...)

reference frame followed by one or more bodies

 
Constrain(...)

velocity_constraints = [...]

u_dependent = [...]

u_auxiliary = [...]

These lists are passed to the KanesMethod object

For more details refer to mechanics/kane and the kane api.
Kane()

KanesMethod(frame, q_ind,u_ind,kd_eqs, q_dependent,configura tion_constraints,u_dep endent,velocity_const raints,acceleration_c onstraints,u_auxiliary)

The KanesMethod object takes a reference frame followed by multiple lists as arguments.

For more details refer to mechanics/kane and the kane api.

Numerical Evaluation and Visualization

Autolev’s CODE Option() command allows one to generate Matlab, C, or Fortran code for numerical evaluation and visualization. Option can be Dynamics, ODE, Nonlinear or Algebraic.

Numerical evaluation for dynamics can be achieved using PyDy. One can pass in the KanesMethod object to the System class along with the values for the constants, specifieds, initial conditions and time steps. The equations of motion can then be integrated. The plotting is achieved using matlplotlib. Here is an example from the PyDy documentation on how it is done:

from numpy import array, linspace, sin
from pydy.system import System

sys = System(kane,
                  constants = {mass: 1.0, stiffness: 1.0,
                               damping: 0.2, gravity: 9.8},
                  specifieds = {force: lambda x, t: sin(t)},
                  initial_conditions = {position: 0.1, speed:-1.0},
                  times = linspace(0.0, 10.0, 1000))

y = sys.integrate()

import matplotlib.pyplot as plt
plt.plot(sys.times, y)
plt.legend((str(position), str(speed)))
plt.show()

For information on all the things PyDy can accomplish refer to the PyDy documentation.

The tools in the PyDy workflow are :

  • SymPy: SymPy is a Python library for
    symbolic computation. It provides computer algebra capabilities either as a standalone application, as a library to other applications, or live on the web as SymPy Live or SymPy Gamma.
  • NumPy: NumPy is a library for the
    Python programming language, adding support for large, multi-dimensional arrays and matrices, along with a large collection of high-level mathematical functions to operate on these arrays.
  • SciPy: SciPy is an open source
    Python library used for scientific computing and technical computing. SciPy contains modules for optimization, linear algebra, integration, interpolation, special functions, FFT, signal and image processing, ODE solvers and other tasks common in science and engineering.
  • IPython: IPython is a command shell
    for interactive computing in multiple programming languages, originally developed for the Python programming language, that offers introspection, rich media, shell syntax, tab completion, and history.
  • Theano: Theano is
    a numerical computation library for Python. In Theano, computations are expressed using a NumPy-esque syntax and compiled to run efficiently on either CPU or GPU architectures
  • Cython: Cython is a superset of the
    Python programming language, designed to give C-like performance with code that is mostly written in Python. Cython is a compiled language that generates CPython extension modules.
  • matplotlib: matplotlib is a
    plotting library for the Python programming language and its numerical mathematics extension NumPy.

One will be able to write code equivalent to the Matlab, C or Fortran code generated by Autolev using these scientific computing tools. It is recommended to go over these modules to gain an understanding of scientific computing with Python.

Links

SymPy tutorial

SymPy documentation

SymPy physics vector documentation

SymPy mechanics documentation

PyDy documentation