Upgrade to Pro — share decks privately, control downloads, hide ads and more …

Python for the Dynamicist

Jim Crist
March 11, 2015

Python for the Dynamicist

A presentation given for ME 8287: Advanced Robotics with Medical Applications at the University of Minnesota on using Python for rigid body dynamics analysis and control.

Materials can be found here: https://github.com/jcrist/talks/tree/master/dynamics_talk

Demo can be downloaded, or you can view it on nbviewer here: http://nbviewer.ipython.org/github/jcrist/talks/blob/master/dynamics_talk/dynamics_talk_demo.ipynb

Jim Crist

March 11, 2015
Tweet

More Decks by Jim Crist

Other Decks in Research

Transcript

  1. Overview • A brief Python Overview • Introduction to SymPy

    • Multibody Dynamics Review • Examples!
  2. Core Data Types • Integers 1, 2, 3 • Floats

    1.1, 2.5, 3.14 • Strings "apple", 'apple' • Constants True, False, None
  3. Core Data Types • Lists ["apple", "orange", "pear"] • Tuples

    ("apple", "orange", "pear") • Dictionaries {"apple": 5, "orange": 6, "pear": 4} • Sets >>> {1, 2, 3, 3} {1, 2, 3}
  4. Functions def count_chars(string): """Return the number of characters in a

    string.""" count = 0 for char in string: count += 1 return count
  5. Why use Python??? • Free • General purpose language •

    Great Standard Library • Fast enough (usually) • Easily extensible • Active community
  6. Active Community • Someone has probably already written the code

    you need • Python Package Index (PyPI) • https://pypi.python.org/pypi • Use pip or conda to install • As of 3/10/15, there were 56245 packages!
  7. Scientific Python Ecosystem • NumPy • Base array manipulation package

    • SciPy • Huge number of scientific routines • Wrappers for many Fortran libraries • Matplotlib • Plotting library • IPython • A better interactive environment • Pandas • Statistical analysis • SymPy • Computer algebra
  8. Numeric Arrays in Python • Not a core datatype, added

    with the NumPy package • Syntax: >>> from numpy import array >>> array([[1, 2], [3, 4]]) array([[1, 2], [3, 4]])
  9. Numeric “Gotchas” • Integer division by default in Python 2

    >>> 4/3 1 >>> 4/3. 1.33333333333333 • In Python 3, automatic promotion to float >>> 4/3 1.33333333333333 >>> 4//3 1
  10. Numeric “Gotchas” • Array multiplication is element-wise >>> a =

    array([[1, 2], [3, 4]]) >>> a * a array([[1, 4], [9, 16]]) • Use dot to do matrix multiplication (dot product) >>> a.dot(a) array([[7, 10], [15, 22]])
  11. To Learn More • Go through the Official Python Tutorial

    • Installation is easy on Mac/Linux • Install with system package manager (brew, apt-get, pacman…) • Can be a bit of a pain on Windows • Recommend Anaconda distribution
  12. What is SymPy??? • A Computer Algebra System (CAS) written

    in pure Python • Allows you to create expressions, and do things with them • SymPy knows about all the math needed
  13. Create Expressions >>> from sympy import * # Create Symbols

    >>> a, b, c = symbols('a, b, c') # Combine Symbols into expressions >>> expr = sin(a) + cos(b**2)*c >>> expr sin(a) + cos(b**2)*c
  14. Manipulate Expressions # subs performs substitutions on expressions >>> expr

    = expr.subs(b, a**(1/2)) >>> expr sin(a) + sin(a)*c # simplify reduces the size of expressions >>> expr = expr.simplify() >>> expr sin(a)*(1 + c)
  15. Solve Expressions # Create another symbol x >>> x =

    symbols('x') # Quadratic equation >>> lhs = a*x**2 + b*x + c >>> expr = Eq(lhs, 0) >>> expr a*x**2 + b*x + c == 0 # Solve returns solutions to expressions >>> solve(expr, x) [(-b + sqrt(-4*a*c + b**2))/(2*a), -(b + sqrt(-4*a*c + b**2))/(2*a)]
  16. Submodules (just a few…) • Calculus • Linear Algebra •

    Set Theory • Combinatorics • Differential Geometry • Tensors • Classical Mechanics • etc…
  17. Submodules (just a few…) • Calculus • Linear Algebra •

    Set Theory • Combinatorics • Differential Geometry • Tensors • Classical Mechanics • etc… What this talk is about
  18. Classical Mechanics • Describes the motion of macroscopic objects •

    Concerned with solving = = −1 ∙ • Many ways to do so: • Newton • Lagrange • Kane
  19. Constraints • If coordinates and speeds are nonminimal, need constraints

    • Holonomic • Relations between coordinates • 1 , 2 , 3 , … , , = 0 • Nonholonomic • “Not-holonomic” (i.e. cannot be put in form above)
  20. Constraints • Pendulum has 1 DOF, but 2 coordinates •

    Length is always • 1 , 2 , = 1 2 + 2 2 − 2 = 0
  21. Formulations • Lagrange • Enforces constraints with Lagrange Multipliers ()

    • Energy based formulation • Kane • Avoids need to compute Lagrange Multipliers • Uses generalized forces ( and ∗) • Can be simpler to compute for large systems
  22. Equations of Motion ℎ , = 0×1 ℎ , ,

    = 0×1 , , , = 0×1 , , , , , = 0 −+ ×1 with , ∈ ℝ , ∈ ℝ ∈ ℝ ∈ ℝ
  23. Equations of Motion ℎ , = 0×1 ℎ , ,

    = 0×1 , , , = 0×1 , , , , , = 0 −+ ×1 with , ∈ ℝ , ∈ ℝ ∈ ℝ ∈ ℝ Holonomic Constraints Nonholonomic Constraints Kinematic Equations Dynamic Equations Generalized Coordinates Generalized Speeds System Inputs Lagrange Multipliers
  24. The PyDy Workflow 1. Specify system symbolically 2. Solve for

    EOM using either Kane’s or Lagrange’s method 3. Analyze resulting EOM symbolically 4. Compile equations, simulate system numerically 5. Visualize simulation results
  25. Analysis - Linearization • For systems with constraints, need to

    be careful! • Can’t just take Jacobian, as states are interdependent • i.e. 1 , , instead of 1 ()
  26. Linearization • Routine implemented to robustly linearize all cases •

    Works for formulations from Kane or Lagrange • Handles holonomic and nonholonomic constraints • Applied to Whipple bicycle model [1], matches to 14 digits the benchmark eigenvalues [1] Meijaard, J.P., Papadopoulos, J.M., Ruina, A., Schwab, A.L.: Linearized dynamics equations for the balance and steer of a bicycle: a benchmark and review. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 463(2084), 1955–1982(2007)
  27. Printers Symbolic Expression Fast Numeric Code Generators Wrappers • C

    • Fortran • Matlab/Octave • Javascript • Mathematica • Theano • C • Fortran • Matlab/Octave • Cython • F2Py • Numpy ufuncs Code Generation
  28. Caveats… • C/Fortran use Floating Point numbers • Poorly conditioned

    expressions can give you unrealistic results • No universal way to guarantee well conditioned… • Decent solution: • Simulate at several points with multiple precision numbers (GMP) • Repeat computation with doubles, compare for accuracy • Large systems -> large expressions • Can lead to slow derivation • Most systems are fine (have done 25 link inverted pendulum easily)
  29. Why? • Open Source • Other symbolic solvers (AutoLev*, SD/FAST)

    are proprietary • Fast • Slower generation than numeric methods • Faster simulation • Correct • No math errors for hand-calculated EOM • Automatic translation to C/Fortran -> no typing errors *Now defunct
  30. Future Work • O(n) articulated methods • Featherstone/Jain • High-level

    specification • Specify link/joint parameters to easily build up model • Finish work on LAPACK/BLAS/ODEPACK integration • Will remove Python layer bottleneck for simulation • Integration with Ipopt/CasADi for optimization problems