{
 "metadata": {
  "name": "",
  "signature": "sha256:446bed56e98f63afe6287b9050f148a7b7e5efcf1fc1b8d572afaf8d4af731f4"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import sympy\n",
      "import sympy.physics.mechanics as mech\n",
      "import copy\n",
      "import os\n",
      "import scipy.integrate\n",
      "import collections\n",
      "import matplotlib.pyplot as plt\n",
      "import sympy_utils\n",
      "import numpy as np\n",
      "%matplotlib inline\n",
      "sympy.init_printing()\n",
      "sympy.physics.vector.init_vprinting()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 18
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x, y, z, phi, theta, psi, t, J_x, J_y, J_z, Jxz, L, D, Y, M, \\\n",
      "    alpha, beta, V, h, g, P, Q, R, l, m, n = \\\n",
      "    sympy.symbols('x, y, z, phi, theta, psi, t, J_x, J_y, J_z, Jxz, L, D, Y, M, alpha, beta, V, h, g, P, Q, R, l, m, n')\n",
      "frame_i = mech.ReferenceFrame('i')\n",
      "frame_b = frame_i.orientnew('b', 'Body', (phi(t), theta(t), psi(t)), '321')\n",
      "frame_b2 = copy.copy(frame_b)\n",
      "frame_b.set_ang_vel(frame_i, P(t)*frame_b.x + Q(t)*frame_b.y + R(t)*frame_b.z)\n",
      "frame_w = frame_b.orientnew('w', 'Body', (beta(t), -alpha(t), 0), '321')\n",
      "point_cm = mech.Point('cm')\n",
      "point_cm.set_vel(frame_b, 0)\n",
      "point_cm.set_vel(frame_i, V(t)*frame_w.x)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 19
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "w_ib = frame_i.orientnew('b', 'Body', (phi(t), theta(t), psi(t)), '321').ang_vel_in(frame_i)\n",
      "w_ib"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$(- \\operatorname{sin}\\left(\\theta\\right) \\dot{\\phi} + \\dot{\\psi})\\mathbf{\\hat{b}_x} + (\\operatorname{sin}\\left(\\psi\\right) \\operatorname{cos}\\left(\\theta\\right) \\dot{\\phi} + \\operatorname{cos}\\left(\\psi\\right) \\dot{\\theta})\\mathbf{\\hat{b}_y} + (- \\operatorname{sin}\\left(\\psi\\right) \\dot{\\theta} + \\operatorname{cos}\\left(\\psi\\right) \\operatorname{cos}\\left(\\theta\\right) \\dot{\\phi})\\mathbf{\\hat{b}_z}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAukAAAAcBAMAAAAuKjy7AAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMmYiu80QdonvRN2Z\nVKvu110NAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAJ/UlEQVRoBd1Zf4hcVxX+3ryd2dnZ3exW0LBq\nu6NR24jgGAsSauiStgZBzIIZrSayo1WMNLqDoglU7dMSG2ptt4ZGLa0ZW/+wFMlqqVisOlpJKrUy\nYKFqC7sVURTJpklssmuy4/lxf703b9NQFhpz/3j33Hu+853vnLlz32QDrMoorQrLapFcWGp8Vaut\n66oPeu5X3rqw1Ph+vExdg/Oegiy3rNz9aC3lWaWFS3BefBa9Cmos1XnlfSnQLw3g5eqankxlcMvp\nZnkm5VmlhUtwXnwWvQpqLNV55X0JUOk/BrAKulKp6MJqpDZe0cWFpca34kLV5RVe7FZUXblC9pXe\nCLzBQjK396i4ou4p6w/nV4ULsaOqTBolZuYRVSVfgRAyYjPnTspPUAv/eJWFrqAmquaS6KajsoUS\n1YojEl9EfovOQJmNCKwsUBVsj3ebHnilN8Uay/h+8Jk21lCQjMTMZqIlu8byuh5PCChLR5smStzp\nB2vhfM+b7XLanVoZ/oQ2FV55z3Mrq0G2zOjbns1TmUKZasWhTCzNtSWFZTYhUFkAQ8mODgadzf7m\nG088B/miA3gRmDR7gY93eEmuYl7Xiw1y0iecyCQPk4p3LKF3smXzFc32ubpu+BOCKnysMdJYUQ1T\np0Z02i8DKtUlVN6ftpRJpOVWwWxCEFbB9lTTEfUnzuwx2DdYxVlgj/pKGneDQcZVMsiV6bq6TYiB\n8mRS2ajAY0yXzyDxKXGUar1QK0nIFP4TjHRy1Gi0ZcxjMtUFuoQKuCEPbZg2sK+3RrMpBDYnV8F2\n0PWB3JIIR4N964HjwKd1ozIq8zO6QoFjyZXpurr3G5CfTCob5R3WcvlKM7r1fZmGmxYQzsovZAIv\nzGCkmqNGo89VJgIqKVSpAFtnmFaaQhvreNO0JeVnNiUIq2A76PpD6RDgHaNuh32PIaKzXuzIpvnG\nWzWyJFfx1N76YRelYiuzupGho00b5QOs5fN9VLcSmXK7bvhVEcMH2xgnKT1qNLqnzPgSTUHPkEoK\nVaoVum6YEg43bWHTDWYzBGEVZE/tqe9oKO52Bydj7TNfiLfORx/+9e5/8Db5KkvoX6JOtXl9ydCm\nGs/S9Z1HqH83iavYrUaL4mGvujmilw5BFLBr3zdRunVPi4Ec6POJqmj08T/ytvYtA2Z+S8aBONTC\nzS0WmlGj0ULIZDyiL++evHIZr7537VeatAyp2DZUppDSE/s62Ln7NokiJzGVbm8i6afuCjpHmGoR\nWTBVUNjUcYzr1wqbiciO0gQG0TeP4hZczU0k35rlo/88Rl+ZSVoOdTbgA4yVridz2ISbG+wqngHu\nZIQOcfdVkUOHIArxfRirva6GLQKkWJ/vzcz0VFwtsgzpWxZM/I4MDF84evSgCM2o0a6HZeITwCze\nD2xrxicoMqSSQpXK1LmxE+2PHkZfU6JEZHn4WNzpf1HbkifMEIRVkD3VRtn8O/UXlNaOygE6vtT1\nPrkiAfINTmKE+jk8S5jNdJndwVhua9x5Em/BgLi4zjnWr0O6PthBDl0YheIk5kZ/Rj9zXst5afh8\nz9Kqr1qolefJkL5lwcTvyMDwRwAqiYRm1GjXwzLxYAMJd/298s5CSCWFKpXpeh3DLwwmiF6QKBG5\naWi20MCMtiVHmNEismwVJJG6PtzlYwTcQ1+Y7T+n8cgo8KstNek6tbqpvpGEb0tEx+g1PEtan+QY\nbmsJ36HlkLi4zvGz7JEhXR9o5dGFUZhuUzSFjSSclwZpMfnGaXUTfRRDs2RI37Jg4ndkYPhp4SKh\nGTXa9VSZgyc/D+76dcApigypuFBDZeqk2qQFpyVKRDYGEjolbW1LjjBLEFZBNnU96naITkhklsfO\nqTPS9QnX9ekOfkcuFkOviATbsLdef7Re34HKCVrSj6JMndbNXe+lC6NwqEkdXaJOTzCQBrXG5GO9\nj1HXixPlen373fV6OwtmfiuBu146i6GZHjU2mqn9KO3ttrnr17quOyouVKlMIRX+XOZawEmJ0oYd\nqpXBv6IZnSPMaJHDoFXQB5fqevjVi2t4e4tvGNN18k235Djwi3mkGldjOdFymMuzcQdfl18AfLoO\n9d4wOXQIovxZf5zzUnE+H/9QOI4NuLpDRvqsGzBdC44MBC8tyUkhoRk1OTfMO1FYCroeUnGhhipz\n1s9IlIjEN+iVVm7orx931r0wSxBWQTad9Yq5YcLXzHCCNdWg6+QbaJVJh7xNi81Cra9NC32bDk0U\nGsOt4P31bLdx8KfOTa+oHDoEUXyvo3MvfWA/4rwU6fPxe+hD9Bp5H/NJ37Jg4ndk8ja9h3oRqFm/\n2BiXA65dD8vEn4C/Bl0PqeRtqlSmzjoJGKB+LUmUiMRv6ZX2F83WUwWzGYKwCrKp6/Q27W6jW0Ze\nj0TAY/gA+mp00O1ZJ99Qc2ONPENt/mDL2gXtemWmH29VV3FxNPpuDVsf0O+xfhXm8+gQRCH+Fsq1\nsVFcL3mJyuf7Pa0eQqIfsvQtC+b3rCUDw4/E++lJQlVN/F8MNGhDP7NUmfhzA9fw9WJumJCKCzVU\npusbq7glugvFRKJEJB5Es9AiIKPzhKkWkWWrIIlTd2Fsf+lE+STX5kfla5+7LVpY/tjC8qVbv0fb\n5CsdOcz+Yoser3/TFaM0GTV47unL+QMptlA8/IfrErr4Fzu0Nm76qubQIYgCPnvkFsRP3NFiIAf6\nfB+hVfSbzV/lXdO3DJivAksGhu/6UoOeRadmofU2WtvosEz88PIrWpd1/35Z9/mnF+f5OBHMV2Op\nTCGlf++p4d37/iVRBGSm4ac20/tYsvVUwWyqRWTZKkjijl1/ux6l9twBalWNw/OH920UgPztgSw5\nzMClsqkuNuPFhuyo+2Gxw4fSZaMChMsXz/BuXFWf3hFq+6fwK5nCxeXVjE+qAo121J7AWwGVjye3\nqdMD2RKmqCObKbSFucLDKrzEqPuaBgpti+6dve9Wca4zkPt1/qJM6mKzb0GudahbT5oi9al02agA\n4fKVm7xbMAeiUg0wzhR+JVO4eAI1y/Oyo9GO2sUHRkDl48lv6gyQZApTf0M2U2gLc4WHVXiJN55l\n+xqLzpmtrzQpziQNkaVxsedohW4sN/pMz9yGSZXwRhAV+J2WT8qm/o0l9Ie28Ceyo3A2A97I/FwQ\nxLnL9FRBvInrnbgpKi0f7QoPq/AStyY3dtDz5/4wjf2vAENEn1E4ZOlyYP0pdO/z/njC28YSukxU\nCmTz6afdn/JlF8KviuzhoK9bzcOu9SZZljq1aRaeKozPQ/IeM6m0fLQrPKzCSax0t3db/vZkwuyw\nN+u6rMOvV3bZt4DHWjrzd1LvsJYB2KvFbufPjj+AezXvwo9TYT53alsXjsrH56DMlmdaAW3YAlnu\nslyZ9eLw3Ll2/uIo5P+qiqt2XMBy/wcsZHK/MjujRgAAAABJRU5ErkJggg==\n",
       "prompt_number": 20,
       "text": [
        "(-sin(\u03b8)\u22c5\u03c6\u0307 + \u03c8\u0307) \u001b[94m\u001b[1mb_x\u001b[0;0m\u001b[0;0m + (sin(\u03c8)\u22c5cos(\u03b8)\u22c5\u03c6\u0307 + cos(\u03c8)\u22c5\u03b8\u0307) \u001b[94m\u001b[1mb_y\u001b[0;0m\u001b[0;0m + (-sin(\u03c8)\u22c5\u03b8\u0307 + cos(\u03c8)\u22c5cos(\u03b8)\u22c5\u03c6\u0307) \u001b[94m\u001b[1mb_z\u001b[0;0m\u001b[0;0m"
       ]
      }
     ],
     "prompt_number": 20
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "P_eq = w_ib.dot(frame_b.x).simplify()\n",
      "Q_eq = w_ib.dot(frame_b.y).simplify()\n",
      "R_eq = w_ib.dot(frame_b.z).simplify()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 21
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "euler_sol = sympy.solve([P(t) - P_eq, Q(t) - Q_eq, R(t) - R_eq],\n",
      "                        [phi(t).diff(t), theta(t).diff(t), psi(t).diff(t)])\n",
      "euler_sol"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{Bmatrix}\\dot{\\phi} : \\frac{1}{\\operatorname{cos}\\left(\\theta\\right)} \\left(Q \\operatorname{sin}\\left(\\psi\\right) + R \\operatorname{cos}\\left(\\psi\\right)\\right), & \\dot{\\psi} : P + Q \\operatorname{sin}\\left(\\psi\\right) \\operatorname{tan}\\left(\\theta\\right) + R \\operatorname{cos}\\left(\\psi\\right) \\operatorname{tan}\\left(\\theta\\right), & \\dot{\\theta} : Q \\operatorname{cos}\\left(\\psi\\right) - R \\operatorname{sin}\\left(\\psi\\right)\\end{Bmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAA8YAAAAmBAMAAAD0E5/tAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAVLsyIs3dRBCZ76tm\niXbIwtSaAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAM3ElEQVR4Ae1abYxcZRV+5mtnZndmd4KKWiMd\npP5QsB0kgHyEvSB/qCVdq62g0Y6Gb9SOmPADPzqoEClCR0QiEtqlQGxTNSuoiAU6UYxALY4JCE1p\nO8pHSKDb0iJQFro+5/26987cmVlIadqyJ7nvfd9zzvuc5z3nfrx3dgEgdhSbaTkUM7CholeVXGQ6\nh+Ii3+VrGlxbVRm4euRdnohDefnLdHFX6lIfyit9964t/YZa+/k6A31jkZn4ZKR2fys7kNvfNPZV\nvEv3FVBvnD3KZZV2XB75yD7pj71h9oNHNLn9EPgdCZH53zsCGwm6VmlNjSM9gMM76KfVB0cGdHWn\na3xwVOvtsZyu8dvL28E0K1xjL5J627P61Eg3rUzU5ZwAciWtmEqrJ4Ez/kp3OSLkLeC1ztaEDXyr\nUcZCIPOPjqGjpmCqoJGJkGBTkZ5BIkDcHBM5XONGxIz293GySLfkpi33We+Tm7YHnKG6ebZ3+0rT\n2zC5bt0xbVoq9CQw/AALIUeEhMi9MHn0XedWW73OvnVz3egSv/WNijAsvK/3e0LgsitrEaE7cp46\naHsidDA/PALEA1pmuaiGJjEhE3LrJ9bdVde66NWaErydGqeqBH68hIERHQDDgfSfo3RS45SxBk67\ngPnFwNh29SQpAgTT4lqzPgeCUMF9aXx12AH5zwF/MLrEm75REVY1NtAV36Z7JJBYjdeiQnfi3Aqa\nbsWEBY1IhA7mZgSJO6V0WoOEjENFJMlOpMNqdeRQjTNNNaG1aX1Wb6fDVUU+3Na0egKDuhLzaDFd\n4EfWLbGXN+luO/LPxjNZp0rA5WgTZXVaWVSMcCIWP7OAg2U1pQs1Gs/BI9PqIwT66njJD20x0Ymz\nIemDXhMKyYEDdYnwPXQwO+5IvC2InaHOw01gfSmkkkFgtTryKuWiW6TbJ4h5hvLxmzvYXSTD9b7O\n9vo91XuQbWZUdYEjzBmDrG/fy3bkn82knEfVMnP4VtNTVqcdfJVJ1B/3Dj9Vpnlgp/NxHSHMDYLH\nRuCRq0kbECHwCPCGsdPSkzNaQR8O4KmuA3WJ8D10MDvuSLwtiJ2hzjtYriVeSCWDADEdWVd3Qjvm\n2/xFcdaSB0L6dBnIM8XANnNNfNqcqbpc9PqxiM/qvp+v+CgwszW/dDKTVPhUgU8oHm0SJhffyXLy\nvhOx9dhW4SBvnhPJ45RNGiFMcfDtHITAxUgQz4a2mOjEuRU0Zi44FUo1PqhNhG/Twew4RNwq5dwa\nJGgDPs5HmYnaYbW6BOcrKJ2s47KnezLs8TNbvsa36og4rmRSTzvi2eScscQXfzauroRNvH83NdEY\nLAPsK3H5GhpD7JtUZW7dUsB1408g8Z9xARLH61awCM+xEkSXIyzWyufAmLIMcd5842bwk+qnoyyv\nPoE9YwLHrznthSa9Bc4CSD+3/tLfITe+FUHW6VchzwYbuifnVtANr6wr4NHtTRdVVmVAbSLcmoze\npNonjmu2/AuZrdsrktYI5g5AdXhnHt7sulpdgm1VuvcxYUC2MA+f5zlzp4w6S7zOO5gHMKeETBF9\niI8hNRvzPaq4bc7HdiYLg9y//E18KC5fww/8l49DYG4hcUfiPMSbPwDKHMteuzETp2NHFTlykSMs\n1sqHqbZxm5f+hvEx+HmFHed9rGCPBhY0k/JmEMIWQEH31ZgdLxtiPTCx7s87/dA9ObeBXsy76jzc\n4qLKqgyoTYRbk9GbVPvEk7fhZO8ED7NVWiOYOwR2EhPjc7iybqvVJVg2xg36gipdyelBbJZOD+kr\n2NfA67xSV/P2Y435SBuqc+JXgdOz5VwVo/43l8sXS/iIR6eFiO3uayCx+4qq2ktzUrJwEx5CP/NU\n1gfdfHFWp9px/6YVnhkZfP2mHxoFFCxrzIeZ1J2EHYDAQ2o8o5R+Lcia17k8G5SdLj05t4FKjW/H\n0qqNKqkwoIcTMCRGb3Q+8dQIZpZ46Q6fKGmNYB5E4W2Uu4OKLqvVP0MnL61iBovGHUyZ+b7JB5mM\nFNr7K8Ck+Kl97ZLZnqoxU9Sk7ldAtb9BfkzjsDhRXL6Y8yHqM3uoHGbQN/t2PSMOnJTBbxg+y6tz\npz5E78RZnWZO1XUt/hBn8n3fZGoEljX+KLCXKhJ2AAKvaozTNr2MuM96qKEoKTtdenJuA2WNkdi+\n0rNRZVUG1CaCHlqM3o4c8eU1YvDFOdSQtEYwzyy6gHJ+SWbytoJksstqTQlyt+PaWTKFF3cD8vUR\nJWcqZUG1XB0ukV5fmc11i/eoGhddjbHUy0M+w9TSHl+48JKFC8+lJzDBd+gY731J/Eyi7Mo8Pllj\nn9lA+mWGv75DjZ2VjlpesR0fv29UdBfwULCs8UdcjR2Aq/FVsxK7EfdZLy/gl5yq7D5mZ86ShSBn\n7tmQWeANl2xUWZUBbaux0Qtfik98aZN3DncFQ0VJq76dQkHE3Qp3N0lmtNtqTY0vZ5L/Sb581Cbr\nSb35sijunGlwn9Q0v971FcxX044Sn7EeZ8uz2mSLDyj8m2/WfNVtl909keCDU7bV7j7ecyZyskGX\nSflysoAnmbeyPqgKiLVaVUIyYMXcc+qrST7PoGD9GgthCyDwvJwH+fGX2P3jAOvlFXXPKztd7H3c\nkXMb6MVo8hto+FNVc2XJqgzo5QQMidEbnU/c3cc3SlpJtJV5EGW4rrfVXVarv1gGb+e0PDdISDVz\nXpwViJJ4FfHqOciXxBivA0srZz+dkKmxBgbqgRrL9ukGvlnfz7Pdath8Sf63NTPyPkasv8Zd54eB\n++gok7LFXDXGZcmmSI6wWKvVyhbYicFPv5rcWjpZKCpYv8ZC2AIo6NRYPvMm9wTftjUmgf5Kvsyp\nNnRPzm2gX0JjeQMzTyzZGvugNhHE12KCmZFPnO9jFNYwvd+StOpUh5ibGeq0o4T8XlS7rVaXYCmR\nkNnLJl3O42MyOWJf/X1gBY5Epij2/Bhf3hfV8hs8DmKrEfeYKpst2bRdgWauwvPNPERsvlI7+cho\nHst9dR2P8QpJNWZUMYsOMik9Ooi/85yt6SP8R3RrtfvqlFyVViz+c9fgqd+LUsEy1SbbQtgCCDyy\nXOoqfkk/FGCdbc71xCR2isXsyLkN9CY0+ovYdqq7j7kqAyqJCK3H6G2qHfHkr5H3eJl+WaVVpzrE\nXHEzzWJekXsThW6r1SXgNojCSwc46S8bS3KO+D6ehcxr+DpwlJjTZTanfOK5ivydKP3Ui08k1k/8\ncP3EsXPWUi+PpdhPj3mGJ3xBGorJV3bxRB3Z2U2u997tHs7acg++9+jGCh3Us+yDhz3qsZ+iQg4s\nqbJxYqx8fI1Rl128q+lMFh+Ju7ZeiAT1AnvK5J9Ombz7sNfprggH4IF1j+EzK55/9ukA68yKXwii\nCs1zT85toAP3e8l7x6//io0qqzKgKhHB9Ri9TbUjjmtXPIbkrZsrklayaAsiFI38fJKPwBvGu65W\nl+BONUOqA8xTrWu2bC4Njj+JFzc2UcbgCPh4JCLlPNXGC7mS6oSaPo9DXlyU5Ki0FHtP6FFEqyaB\nNzhlrjmQ9mToRFvdMNgJ4B8Z+YdJRdiHD87VfU1A+hJeJICpFW3t1EF1IlrW04IXTVynuhvzFhgZ\nBojpyKuUk24fDE0YaAw0jsTV6UaykRhB6txb+Hj8kPJ4WLXxsfeF/PUgV+N5sCqDfFNaynf1qXOr\nJuF55bCVrRzIqaFrtNUNg50A/tGxkaDF9BVhH77dQxMQvQrNcwCz3V1ppg6qE9GynhbUaOLoGaQF\nRoaBOTpysMaN0IT5VfmQiD+z5wPVdBH9hVyZf/VUHnFPTtlJ2cK2ySxq5POdf5pS7ZQamQQVPsMa\nyQH8RBpfwuR8fbi3crIUVqhR3ONJAWhopQw1igA1newhZz2IezxPCVQnomU9GsS10cQR9+jRNQjt\nLRL3qFBzTAmCNW7S5stS3kqjGKhdNllJj2HIi9eB45U5WZRT4pKCGrU0Z3A8qHQ2by0OUUOZhKY0\ncU8fQB0haYZGnQYn3BZlUYQVQNyLsrt/UkAne8SsqYPqRNQjMHxVNHH0DuJD2J4/R21owe2liG6t\nkz7P99R9/LQXKyf48aveVN/Rpnlhz9AoWTfDnBfSdx24Sfw5Fero6v6WjY6wwEeJJdDJHjVnqqBv\nJRFtcXoGaZsR2FeZyKuUj27D7gO1TPMYzL2yhmKmjJQne+qNYZfp0UGRgc41xqbtyG+5MbaV3zkX\nIXFzheu556BY0zTJcAa61Djg+KLpjwR0092DJQO6xl/rQTeu7bF6D79p84GYAV3jRT2opevK4doe\nbtPmAzED6k9A/BU58jMoQPg9qt/ld4iA73T3wMqA+Q+37EsHFq1pNvswAzvqGuyRc/Yh6DTUgZSB\nEy60bN5rO9PnQywDx8l6/g/J9UUqxVSUBQAAAABJRU5ErkJggg==\n",
       "prompt_number": 22,
       "text": [
        "\u23a7    Q\u22c5sin(\u03c8) + R\u22c5cos(\u03c8)                                                      \n",
        "\u23a8\u03c6\u0307: \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500, \u03c8\u0307: P + Q\u22c5sin(\u03c8)\u22c5tan(\u03b8) + R\u22c5cos(\u03c8)\u22c5tan(\u03b8), \u03b8\u0307: Q\u22c5cos\n",
        "\u23a9           cos(\u03b8)                                                            \n",
        "\n",
        "              \u23ab\n",
        "(\u03c8) - R\u22c5sin(\u03c8)\u23ac\n",
        "              \u23ad"
       ]
      }
     ],
     "prompt_number": 22
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "aircraft = mech.RigidBody('aircraft', point_cm, frame_b, M, (mech.inertia(frame_b, J_x, J_y, J_z, 0, 0, Jxz), point_cm))\n",
      "aircraft.set_potential_energy(m*g*h(t))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 23
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "F_A = -D(t)*frame_b.x -Y(t)*frame_b.y -L(t)*frame_b.z\n",
      "M_A = l(t)*frame_b.x + m(t)*frame_b.y + n(t)*frame_b.z"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 24
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "F_T = F_A + M*g*frame_i.z\n",
      "M_T = M_A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 25
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "H_cm_i = aircraft.angular_momentum(point_cm, frame_i)\n",
      "di_H_cm_i = H_cm_i.diff(t, frame_b) + frame_b.ang_vel_in(frame_i).cross(H_cm_i)\n",
      "eom_rot = di_H_cm_i - M_T\n",
      "eom_rot"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$(J_{x} \\dot{P} - J_{y} Q R + Jxz \\dot{R} + \\left(J_{z} R + Jxz P\\right) Q - l)\\mathbf{\\hat{b}_x} + (J_{y} \\dot{Q} + \\left(J_{x} P + Jxz R\\right) R - \\left(J_{z} R + Jxz P\\right) P - m)\\mathbf{\\hat{b}_y} + (J_{y} P Q + J_{z} \\dot{R} + Jxz \\dot{P} - \\left(J_{x} P + Jxz R\\right) Q - n)\\mathbf{\\hat{b}_z}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAABQoAAAAcBAMAAADbxJ9KAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMmYiu80QdonvRN2Z\nVKvu110NAAAACXBIWXMAAA7EAAAOxAGVKw4bAAANU0lEQVR4Ae1abYhcVxl+Z+djJ7s7u4M/lFDs\njh9gglKHEmxNoFlsS+wPu1txg20CuVBrLdZmxNgEq8lYiKkQsxGJYlEcFDVSMVslbdNoM4KJ/VHb\npaL9UTFDm9KCmN1omk1Ss+v7nvO+5+veO3uzaesmeCD3nnOe9+t5zpk7Z+4G4DJvpcu8/ozlX6E0\no4z0l7rZmvVLvcI3pb4rk2b+YJSmTn8nDbmI+TclSIZ8ffsO1Rc2e2phk6VtkY1mBg4LK/F2rRwV\n+4sN30mrefNYGnIR829KkAz5NjfKexc0K72+oMkSN8hEMwOHDEq8XSuH1ZZm4ZsZil76JnheipZ+\nlZdc4RVKM4Ir9Lx7yQv+/wD/CwVytbc+66XkeE9aeQjk5meT0AXSLQAnRTRzadW8w1hcdCctJAda\nNM2LLiTJIVdLms02141YyGoVhSz/4422Dr19/tYTGyPdV1cLwfVP76nh3DXzJ09+zLFwvHk2tPBz\n9EzPnTxBgbzm5NHzXMlgVQ85Ofs8cAoIWJ64C/10kN/90EtuqgCOMVa2e4781vUBLpqrCSnkR8jY\n+IT8dSQD6yH0Tp8+efMxxQNneMSYvWWmaXj4SrmLG6srTXSTPpNWYRRfq5j+FDxkpd90lE9JYjzJ\nF74rA7obqHwbwCs0cxpg7Qh1pBkTmQgsghxDI5BHi6DFgnAl+seSSa69emYAECgm7sIg3c4qDOog\n2jeAIcZYmU1U/fq4aA4UUChGvk/AnyPFQmKwWyLiQW1Ij/TAueJ0JprCI1AKI1llY3VZiDMGYmTT\nKoziaxXqrzL5rHqbarK4l6vInQWoeGsrUOlTaHJPG78HER6kXVuqs5OY8DC0CHMMNwCmq8rYhECp\npQSOIpVso7FJzuDAGAACwfLoaEG6u0bQ/QcqSJ2uAQySh7A76aLb7dLhOxetqgFwKBDOs+JjFEKo\nW8gmwPGGOA/rEbr4LSNN4WGVMtoaZU1dcYhzumLgVEatTAKO4mll9Xd5+ayWqZWBoQ6b9J7H3P9x\n7QUqTuLs4AzCuAP7/42DSoPtxETcAoswx44qwJG6MjYh8FnQEXd9l0ruoaFJzjbLMDUCwS7U0YJ0\nG8hlmi6JMEgesniBLrrNSYfvXLSqBsChQDi/9hIfwx+hbiFbGKhNPKjtaKmR6ruXjDSFh1XKaGuU\nNXXFIc7oioFTgZSSQxk7xEwCjuJpZfVnVN18Vo9o6ECNTQozuNUu8EDdBJpo4ZCevYW9+BFu48Bw\nEROcUy2wCHPcgq680U0IgFgQrqQ4hTFNcs6wtopbcAp34c7xYzyFNx3NT1fGjxX6o30ijHQ4D5lZ\nZfv+RWOncdGqGvwitRTQpm9SGRofwx+nu4QcRXgaSUwpbx6pvnvJRtPwsEoZbY2ypq44xBldMXDK\nl9LkUMaWWGzpXK0c/dNZ7dbQRJ1N6HSytu3aM5THEwPAAK4p7vzK4zQwXMjk7u2rtzWg71ALfy4G\nFmEOfGS82KAATgjcKGEQrqTcBrDJlRfAh/GOQHG+ljuHftx0QX66tRgG4LkWXpJgfR4TxlbZMj7s\nS1tf/vkuXIqHczUALpqqweZQwBFPko9qhj+OuoRERfP4aWRvHnGI/S9ft72xf0srK03h4Sil6WI8\nUlY1U5cLeSyD5fellBw6mCVGCbworlaO/spv/41ILGTFP3Zv15HxwNOEvidkoO4MldXDoYDPwuFX\nt+gHheGCJqVmz+P9Y/CZrc90YhZBjtzc1tEapzAhAGJBuJIe3EY2OfudwDsCRVzECbXN1LyO5qeb\nqBE0WsVLEuwztsrSUedz8MHGQSit++qzeBjmoqkaO6I+tkJN3czxyCiE011CnsUd3lQ80I5OZDTi\n1rypUZptVbCKbDSFh6OU0dYubtLKuSzD5fellBy6REuMls6N4mnl6K/82kgMP60+q98pqIS7V7cd\nr+3+Y10GdBdInQRhCGXZEcFKZSI0yaQclc9U6hDBDegTWAQ5es9ADx+k7ONU5fGDcCWVST6GquRc\n2vfwjgDtwuPyBJJt5qfTB9Bz5Kfr9WEs1mFsld3cAfgGTER/xtf7eexL0VSNHVEfW/+UupGPaoY/\njrqEnHvomS+hBYfkkY6QmzoI+bNArwKy0RQedpmMtrKCsXXRi+uy9MTASrJpRQncKJ5Wjv5ETRHD\nR5nP6mHFumJOQaMRTZQ2HMZ2sIpdgYZQEPWTjs5EQ20ayC4kkxIUR2gu18FLYBHkKOwFwN2jmoTQ\nebwguhKMiIkludSVQx4E0C4cvsCxpCA/3Tyh+le/TubD+JiMjL+zZY7X6G+D9FoA4Gr8J0VTNWZ0\nF6l0+DZY1qJJIB/VDH8c2V0YC3lGW3NIPeKQpdJ5qJyBwZGsNIWHKIWhRVtZwdi6aNFdliyGyAzZ\ntKIEbhRPK9HfIUafLX/xdJ4BFgTgDa2MvQrUj3sH4DD+m8OTYwfK4+Mb9o2Pt3FCmww1sAu/pEtg\nEeTA00keLfD5aUMkBJFKaJFscvLDZxEV4xOx0fx0h8i+fzINdhjvHB8/ND6+kezxIBnRdZYufW28\nSNF6y8iIYGzLWuqmfbCr+WOne8hJ5aU+ZpiAR3oOKjMwMAn0KiAbTeEhSjnaygpKXTHIsgyX35dS\ncoTEdAIbRdRRWon+zIuI0dtbnxU+c/NP0CrplpOHlB460OAMTtFP/Rzu/ONtwtWHzZgMV/EDUSEg\nsKDnurFCeLhmfiLL59XAThCppA9rM8nRnVqhE6lfpvQsPJDwjWzi8TsaenXA9YbVeIz5wbWmqT5u\ngI9c/HP738hXiqZq7Ij62OgbWXxwaPnjoEvINsLYOGRbDeQyMAb9TThQjwjOQBOYh6OUvzwY2Nbl\nQQ5LTwz0CFbOg5mYSO1E8bSaJkZaf+oBEhuqhazw/Nk7q35q7/r99qp6IaSM9UVDfXugP+o7n99V\nXY6rSTtxooGro7cQe1eiHdAb0Yp1QguwOdasX61KKs/ijqWm5DAleEHOKwN9eDfJeW6o1nQOuC/O\nRz96jBAdzaRTdR9oXf/F3LeVXwBDnDEru6OhXrr8tXcG7qcP1qAtumeMQqGqhgIO6deJ9skfvLPl\n8EesS8gWwtg4ZEsN5LKsCZtbcDQ/Zc/xK89Fw+rhFPBQNGXlHKWUGS+PEj555VyW9LrQbUbKdK04\ngRtF1FHEHP1VZCR2PGpaVkOvj16APbR2padwI0Sv9LegaL6alYuGykPRcoAH7oUHX8bZ4gy+KGm8\nG3uKpjaBA1OjcC3kfv3glzGIb2FzwJ9WNAA2ofCzuSlO4OTxg0glA200keTKGr++foZxECieq+a+\nX4fRn6jV4V1oKKm6ex9rl6+pcxa6GTiBMW+ZF+DzHXxGXSjM5Gtw1YoP3GGLVtV4FDBkGa21z8Dz\n+GrH8kesS8gIYWwcMlIDuaxt0duoT5cdmvk3YJky0rvQ8NDLIytnlXKXRwlv63Igj6UEkSpMjnSt\n9Pp7UWSBFTFHfxUVie3sqdrF21yfGFHvJbdtxy2Rg08i402ncXmdpqCPPAJb0ODErifRamDTXA0G\n1pGZVkOZwH1b734Q52HlHXELevepreALAEfnXwL4w1adRIcQ2AsilRRbaMrJtRN+ZJ7GXhH3+7Fn\nb27iN+I5JIBNRzPpdN1Xf/yBFlQT4ATGvGV6d/0T7Uuv3r/tvQDtnn34UJKiiy1EzAj71Og7VftA\nacrlj1iXkOuUMxRbeD86zyPsU0MS2wH+/nWX5nTrQwpLomlWzirlLg8Jn7xyLksTRKXRVeiVS9dK\n424Uo06xRXGs/jQiYqu+5rC6t3Iugv66wuhiz1dmijuvwauq9wnA/5Fjm2whOxP2tIXNkW8FFguH\ngI+yS5AcDICnTiRCLUzHdRem8MOXAMcZ85Yh45Qm1fjwb8xwtelxZ7EhTRybcXhM50mhaTxYKVfb\nQHgXMm4JHbtyl6KV6O8mcFid+iz0tAXLzeSkG96fBH2wurUy5kJ9NXeU1NcWNsdqeKdvtnAIwK84\n1YLkYAA8l02rYyGE6bjuQucvKkIIxxn/mHOl36Qa30I/o3Au38pFPrTYkCaKzViY66jZkIcsj7iw\nUq62gfAuJF5Jd7tyl6KV6O9mMKxKR6awfyNja8uTeAhJbo/2zSjguXn1TEk26jIrOfKvPT/SxSwR\nKsnGD5MbAOBkH/6XJdskHXDdA/P4g8o2gbsxttZ+z0nqAoU6j2741vvd+Sz9lJDi6sC5eUmjQOEh\nNMUDQqXww7EI4b0cl6RVoD8FtqwKZ0v4S3wVF3/fV3bdz93Y7Vdb9Fa46ocxKNOE5OjN7YsyOThG\nZoHD5AaAlbMw71Ym6YDrzh2acgJmYuzau32b1J2F/AgPG8vXe0CGQUpI8XThm2RS3UOaBgyVwkPr\nIoT3cnTbHSZv0DGVB/qTmcHgwOtHzuPHpBY4x4elOv+hNA5lm8mQIzXQ+9KQVMBQSqn7rajm2rQq\nF55P56F8LXwdPOpFEx4pND3bxQ4kx2L8beVx725Y3FrNDDZ/moIs7enLte4UVSfe1UlErjCaiRxx\nsrKimQYt6fnLte4UUddsTAYuf5r/BUXJ2+qypcsmAAAAAElFTkSuQmCC\n",
       "prompt_number": 26,
       "text": [
        "(J\u2093\u22c5P\u0307 - J_y\u22c5Q\u22c5R + Jxz\u22c5R\u0307 + (J_z\u22c5R + Jxz\u22c5P)\u22c5Q - l) \u001b[94m\u001b[1mb_x\u001b[0;0m\u001b[0;0m + (J_y\u22c5Q\u0307 + (J\u2093\u22c5P + Jxz\u22c5R)\u22c5R - (J_z\u22c5R + Jxz\u22c5P)\u22c5P - m) \u001b[94m\u001b[1mb_y\u001b[0;0m\u001b[0;0m + (J_y\u22c5P\u22c5Q + J_z\u22c5R\u0307 + Jxz\u22c5P\u0307 - (J\u2093\u22c5P + Jxz\u22c5R)\u22c5Q - n) \u001b[94m\u001b[1mb_z\u001b[0;0m\u001b[0;0m"
       ]
      }
     ],
     "prompt_number": 26
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "rot_sol = sympy.solve([eom_rot.dot(frame_b.x) - P(t).diff(t),\n",
      "             eom_rot.dot(frame_b.y) - Q(t).diff(t),\n",
      "             eom_rot.dot(frame_b.z) - R(t).diff(t)], [P(t).diff(t), Q(t).diff(t), R(t).diff(t)])\n",
      "rot_sol"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{Bmatrix}\\dot{P} : \\frac{1}{Jxz^{2} - \\left(J_{x} - 1\\right) \\left(J_{z} - 1\\right)} \\left(Jxz \\left(J_{x} P Q - J_{y} P Q + Jxz Q R + n\\right) - \\left(J_{z} - 1\\right) \\left(J_{y} Q R - J_{z} Q R - Jxz P Q + l\\right)\\right), & \\dot{Q} : \\frac{1}{J_{y} - 1} \\left(- J_{x} P R + J_{z} P R + Jxz P^{2} - Jxz R^{2} + m\\right), & \\dot{R} : \\frac{1}{Jxz^{2} - \\left(J_{x} - 1\\right) \\left(J_{z} - 1\\right)} \\left(Jxz \\left(J_{y} Q R - J_{z} Q R - Jxz P Q + l\\right) - \\left(J_{x} - 1\\right) \\left(J_{x} P Q - J_{y} P Q + Jxz Q R + n\\right)\\right)\\end{Bmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAACUIAAAAmBAMAAAD5Um4SAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAVLsyIs3dRBCZ76tm\niXbIwtSaAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAWfUlEQVR4Ae1dD4ycxXV/++92z3e7d2lT1BIa\nrkCUBFp7kwraUlRvURTaJuQuf0xJE5VrG0OAEJ9QazWixF8jTFUj8CUtURFK2EIQjVDJlfypVGN7\n1bSVYmq4ViUqMpy3pCWAY/uAlD81tfvezLyZN/PNt3vr+nzn+BvJ8828efPm9977zez3fbscAACF\nc7HKSx6BPAJ5BFZTBB6d0WjKV5jGagKXY8kjkEfgNI9A48FEReDWqdM8ELn7eQTyCKzGCGzTR9M9\nyWoEl2PKI5BH4DSPQPUNFYCNJgzXDx6O0sTgc/IZJzUCGVl990kFsQKLZfjdC0nO5l7RWZGx19Wq\n9+q1a/89MIbS3dMDz8knnNQIZGT1l74VoviFULD8/foX/mKpi2SrNuImMvyOK2tpzuZe0VmZsQfV\nsuaEOh4IhenjmZXPWfkInBFAWIn7h3+DbwYoMrvZqtW5zEmDDuRsHjRiy66vzyZdH9dieU6PK2yr\nYFJ4QlWSkw/qE7C+tcRVe6g+tUQT/dVyNveP0UnWOLknVKkTd++SuPgkS0ebqQWzABvFb5vrKsBP\nUBhOyo2oIDyhvhzVWl7hAzA5s8QVeqgOJUu00VdtoBMqK+tZ8r6rnziFCBN6o2Lu9yH8iUPY2xLh\nN5DkCXXg2AUPfzTpPTU1OlBOL6Xp9W+9MuubKU9Qf9+ma7W4sfulHW+j9xPlJ59+RItkPbAeTzbm\nRncf2fFwh4Xy+nXZUe0MwGrsSoARfaYF+OHyL+3roMqjx3bs+FmlKivnfm881k0z2WYnEpVnX1ZQ\nDBy5WlYQUSc4ocqvyXkOJa8bD5vV46T1UeM12CqcnbCIr6Hj1mJa1UxpBHwisZ3FZpd0HYTNYdZ5\ngVC+LGzgcMf2SGNTM82EEBWjpSvdPRvuB4S3aTJOyFlyJy9Vz8zvzXxUUvg1pHvVHF0Dvigv/pUP\nom9vkJwCbmos9ZcDq+b5Yr/e7wBjUwDnJwBPNGEEm2EZWM8YYHNjE1B+KTRK/UpKmAGY9C6lHa3R\nBfjrHwb9egUXmZwgXa9Y9/vgsW6a2ZwdniaMji5qKJFgZQUxdULVZ4VBmSReNx426w0nrY8aL8JW\nP80Cdw0dZ4sRVT2pjO6nCs9KDfQSDMLmIOvWbCBfJjZwuCNsgMeZmBaTPoWwG8aWNMrrEsv9kPAm\nTdYJ0nfFJp+Pjr56ei6Dzs4RMVlvR3lClf4X//OXo279JbUGyWljXpmstAPLB3X/t1i8HvXO7sLW\nCYDaV1norr30aq2InhZZc+u7ALubSniL08aWwedkWYBR43N30wmlkfv4ax/CgW2zACWM5Uh4GmPg\n2ziOJYZHj+jauqm7nB07DcXs7fCUhmKASCtSXcpTJ9SacW+40jZdXhfiYbN6JmlSjeGhJauGbSps\nNTgX1VjgOFu0qn7SaMqn1Dy/EnD9gV69QdjsZ91Z9eXLxAaI7BEb7juZmA6UYSqEsUWN2jt3J8z9\nkPAmTdYJYZCalTZWVGJ6Ik2VNinZYjkpcmTBay0KogYjT6jGqyjVPz+wxvo2BsnpmpYyNzYXWP2K\n7h9h8eEZgMOzcAX1d7PQXXvpFboRPS2y5g43ATa1lPCtThtbtbbXBcgCTGoNOqG2qQk+/so0CkcW\nUQNPp6EfKg1ZsfsxPFLPuqmFnB07DcXs7ZquhqLhSCOZQUydUOsTbx6jBF4X4mGzeiZpUo3hoWGr\nZhZhq9/xFtWdwHG2aFX9pNGcL0SsCLiR0QzRIGz2s+4M+vJlYgNE9ogN90NMTAcKfFRiAJuHE+Z+\nSHiTJuuEP89lNaYn0hRk31JY5MiC10sQk/V2lCdUEXfVyIsBiH7dQrufhhvfrJtbOk5Ereq06ptf\nj2J7Hf7bPV7HAzN2W9pTTzjq9MgOOHPn49Y2B7GIIul8kCpRMgArDXVCVcaxHeDfP4Myuv8ttvFW\ncFYpy8q4H8Uj9AL4YLLjpqEuezvZxM8zhKLgCBvSa09MnTN9yR6/a5Nk1sXn7mjYrJ5OmqfG8NCy\nVTOrGKvV+SpC90vouFnYqQZJw9kv+BZUT8CNjGaIBmBzkHVrMJAvExtie4TDXUbuhUwIUFm0qkEn\nlOZ+SHiTJuuEP89lNaYn0uRn31FY5IjBmxUUfrUd9QmlP7Xo0XZyFlWG5gIk2d3yF490s0eDkSd1\nf38Laoe+8dkFgM0PlDq4i2hNvNjbDXz4Lb8OkwgH4J4ZqmXpqSccJb07Dlx2sAvV62YuEebQ2zO6\n2qKJYvVf33NoF0oMQLuaA2xFtqFOKAXdx19G9ADDeLziJ0ch8u4E3afi3BN41IiuAvjq5Rxmx01D\nNfYWXzyoKNoHIWtHqAdZ/fVN5LEr+Gggi0GpXwoqVgiYJmykb/V00gCEGsMzaiIZ7M2Fb/qxhKzI\nQo57/NAWnSonbR8mrbrwdpy7JZEGdNviCPxOazrJIGwWWffQCjkaXi420IufcI9wuEcX7Z6yrvmo\nrFg31AmluB4Sfgz34OSscyKYaJMf0/NJImLkOGlz5Kh883kfONA1+BWajbRmVd844aNt9Qbqb1Nn\nA7VObDFfbOGz8O3wT92NUFv7/N24QrGjlqm01UU/1G6eh/0d6q9rGqm9VNqmSQ+/oR5nCVUqbajN\nj356aAo+d+iuOWeudOTQug6OUzFRHC6uhe0tgH9UQlc5wE7GLXVCjVKkfPz1N0ijiJ9j63c9o9o8\nw1zNqwDrnsQjVCs+fHrxQNmx00iVvX0Y2wRFwaEBW4R6n6yeY+eohn1hwayQMAX5WM8kA6Qaw0OD\nqCaTwd587dgxf1nsVdo+P4xFp2pWr1PSvq9uBdenSCJw9PEbVzyuIrLeg83LxYZe3B+ZSDNBoE17\nq04oxf2Q8Dr5zolgLiffpNPX80kiYmQ5GePK7DndGn5EKSYrSPsTXBQ3MZXDO5/8cku1lqn6HWW3\nhifnv8P+5M1Qg/IciobGlXw7tak0jjx913/yq6L/0TJR99QTmwL16kn9tUILEnifMIdHy+hXjD0T\nxXdhBCabqVcz4ACL9U1TnVDqvYWPX795Gmurx/s9rdREcp/KphbV6J7EQyJTAvicHTuN1NhbfPEA\nBEXBoQFbPHUrjTU+5gkZpWOFhOnIZ/VM0jxvGB6eTui0TAZ7461pOpRfyQ+5sFKRSZunL65iP6pK\nzYot9f+RiaxLtAGbl4sNsT3C4V7TTTNBoE07rU6oM0geEl4fCc4Jf65Nvkmnr+eTRMTIclLmyIAv\njW8EekpVTFaQtiEdyh9K1Mrr9MWiOHYiyyKafUCZLuCNRQL0hRfAL1O1ZoZqOLujLnj/Qe+gAdSn\nq/pusXbF1Vg2NvV4th6Oc5awiXo1qEzQpBJ6ac3R66HXSYrFRDHZPqOYvp5kW2mxqz9MTQeYerjR\nBBAMLxomt3z8YySiryLp/cnYLHWERURo7qvSeKT1FHww2bHTyDCnFVOqoCg4NGCLp07SsUhSFeJr\nPfcYJd7EJtpYJGw4YPVM0tTLN46uSwapyWRYq8p2ynGPH3JhpS6SdnYCdP84Oa4GZBXOivkdCcUS\nRbiSyHoPNi8XG1J7BBFxuOmOUjFBxFWg1VGSnFQnlOJ+SHidfOdEnMpMEqGHi7gTirIvYmQ5KXNk\nwNdqrwI9pSr8ClL5+gTONPl9RWPPqEtPn2dHyi3bFA1l5j1CYJtmAe3/sDqAjtJgdZZqE717Eupg\nGZtWl+uoHtJtJTBVtl59w4Yr7t+wQRnFN1gJThjr0qzP4D9rbmwOykdI+sSGDddt2PBRap6tf8Gk\nQkJ9UyRglvE1dUIZXENt0rga/+Eik3PU8Yp2P4bHUwvhg8mOdQOct402TqWcpk8op+5bT/c+5okY\nJfC69FYtFTZ836Y/TWzSnJqDh4a1mkuGteotqjs6jo4fdmE17CftiwB0/zg5k7Ljz0oNL43NPXju\ns9ah9eWwXGxI7RER7hfQ25AJ4R7zAhKeUDapJk3OCW+aS35aT6SJs29jZDlpcyTAFxZhZEqeUDB6\nH9ymj54Sf/Q5GOXv4qseU4qwmZv4GgnvvMQYdt8NcAleavNZIziI95DlG/ShU3oZP1Php1Con/Le\nO0+buv78r3wd39/MKvFuqg83VdtWWu/mOfhkXK/Q1apaD3WaeHoXyKA1t75jv8qz5zy+Jn4D8TgP\ntRUH+L1XXqZFtlYnVHUa+3QH7XCNLKKIfmhQwk+O8Ks84z5Zi+Jh62n4nB07jVSNt8W5RH+jqOCw\nDXX11M3I5ZvmpQ5mjoo4oQzKhT870AReFyJhE95w0jw1A49zLpJhreq1Xc2OC37IhZXiW406Ju0G\nqP4Q87a+5SyYVmoWgPQ7ZLMcQwuGzdk8l1nHVEfZTCxdJjZwuL30MvcvwIiETHAc1XvMi5c6oRT3\nHeFl8q0TcpbOanUfDCVMEl/PpImz72JkQcscMfjhKRjrJBq/3o6bcfc+PkNL028awlLEUVMqsD3h\n9juoIcYA6G951JsoJJ3oCMrxPVzjqPqK8icbi/AcnRwjpN7Bk6hLXzI/NzJ1DbYVHNgyc/n3Svfh\nuCxab/61Ujuux45qvUJyGBoJHYRzzhweevWjkCirTPa/hsJr9eibcgP483u7aoKr1Ak1OoUCh59w\nVV8tLzQvxlDQKbW/i2wRRViz7kk8rJqGz9mx00jVeDvWmc96U+6ps/V3coOu/FdY6GnJFI2ykXxz\naMaxQsI0YRPecNLoI8VG18AzajIZ7A2vaK/acZD8kBaVnksaPhIMT/1q9Lu81CycKvxOsVmMOTZn\n81xm3UN7ctjA4fbSy9x/A5k8SsQUxXFU7zExhBspwS6/KY8k31JaztJ69bHkYnt0+Ho+SURGLWiZ\nIwa/Zh7/Q6h58aa8QUdAXd2qV/iGXXyXNzYuQJ3B7WqXWt4Y/VdetQmAz2aO4MA+2lO16/GEf7G4\nWO7ARXsf+00U1+foofVWrJuVDvavSrBCtz85W3+0pZquUnrlP5luTMT12FFtb8v4Ong/lD711H/h\nVmNzV2H+jpa0YyaK5VehPvVtgLvcQqplAQP8cTCkf7E5PItii1/jevZP4alvoLiyiGd/913ePO2+\nshbFw8pp+JwdO41Ujbdjf9jFW2mEouBsn2IreBXq26xcf3XLWiat9GxgikZZgo9jn9eFSNiEN5w0\nT83AM2oyGdYqL8lX7bjHD7mwUnNJwzcWlfl5gD0k9/wWODL8DticigmyudpFs3Gei6xnsXk52cDh\nFum1bIBX7zBMIPSmWI6aPcZydVUnlOK+JbyffKa0nKWz+mub4RlHEk/PpEnryYxa0DKzvHEnZ+CJ\n0aZmsoK0pYOr1o5iNXzVS128UHG/IJlsKYGqah/hdrFFLTmmvwY7F+C8zBEc2Iz/Dh4YxwV3PXfw\nnwFmR++fQVF1GnfSwt9gC8ZaAHceW0tNfI3+c8/OQFO3udZ6I/Mj3bgeO6r1bjt0x1MdGIY9dBAa\nc3ceewR/hHxIGzRRRATltYiL1GSxgKFMQGUZ/rtX0E7Fx69wlR5euAZKGM8juPTarpxk3NfWYnhY\nOQXfZUdGxXjb+BLOIygKzpD9pCFrTt1mtYzf5opituqWxMlUkgBQz64bC5tOpvKGk+apcTK0OZEM\na9UtaVracckPz6LScklr4C/OdjbNb8o9v8WsuN8BmyMxOReKLVovynPJWolWypePDb25//y8YQKh\nNyXYYyxW17u/dn3HcN8R3ku+obQ3yyR/J+xyJPH0TJoiOz7GfOYKIrj0Pwx+tR3PV6t+1V/b9fDD\niR5IyxtvmYHbEpYPqYb64GKR+i/lcdtPkyA6gvKhFo2myyesaL1/IhXH8TiNlKGZNZ7c6bGjkUl4\nY54yx1HU6uV2MM0Bvgx+MRhT3d9QtcXPuH6mzx9CMdbSeGJrpGRumvSWoGg4v+/PcOpGPtr2FMwJ\nZd6liqHSYkn0ZNMPW0Zs3Ie6nHoC2v7qZPC3ldXA79RK5HcWmyMxOQQ9eI7GbdaDhax8pdng4bKo\nINhjRktzf6jFk8LkZ1H6GriPp6ir00unySmmOMmPA6xCTNaQ7leiB3kgvO7A1zj4NmL4rAWo/fnt\nPHqhauzgrroSz98CpSnqxEZIPjpLdbp8x4jKi9+FeTlcnPsJ2bXtoeYW26aG06t2vAG/49RY/kfc\nUNd61+sKwOWdZ00EY6q7oGrGD4zrgoKKQ2wGydhaGk/WDE/upklvCYqG0/G0RXCMfGTOUzAn1LAX\nedSYrE/XPUXX8cLG3rhhbkl4LDsBV291slfTd42dPrbR70w2R2LyFsjmOa1ksx4sa+UrzQYPF6NK\n7TGjpblvd2gq+VmU/r3qoreO00ulSeg5CrMw4AoxWUO6V2nompXF9SF8XPk49mvjMHTsFR64UTUe\n4q66Es9vhOoEdWIjJNfPgKrlVcWW7tZ23Xwo0U1dDx+jX/qkS3XX2zxhlp6nhI8sGeZYbSs37PU8\n02qU7k+s0DVq+iAqtoyIcd1zrOmU0i221g9PeqaSRKcRFA2n1vHnpdTXdD0Fc0KV2p4U4LbvLzwX\niKJd9iY6eFKEo3O0TOh3amn0O5PNkZjciHymksXmYotG06XYMrIVZoOPjFGl9phRM9xnwqeSn0Xp\nm57Re4BXy9LjcX1NcdIf1kzWkHqfUDU6IOithfmqvXAOlnfATWRPjeGdvBapp7wLQf3J6OgITYFL\nVZ2qyhMpkRaUrhuPjtxCXwuKkqUnVKjZT43T46ZlAGaFYku1LH7GddFfskbPaz88GZOj04otvJXE\nf3jjR5UoKfXtTTUqM0cC9SJQzDt1mkPKodDvFH7ldwabIzG5sAfPybbNerCQla8wG3xYFpUvtj3D\n/UzCZ1C61ip2rA1qZOh5Ov13YrGFEzSkXidUqT3aRmP4NqI8U0pwChf8tYga+8DCoyxSJ9RNUJrQ\nI8PXfiQYUd1yx0q9xvu9Xv/Ozqp/cPefsRSN0VZKKwuwUfwHc2X8y4MrhSomICgMJzZuZRfja8LC\nAkxYgbmHgq1Wcqo19i4FMPmdxWYaKx/8wZy1QzG5iX4VpXiewWbOup1lGixfYTYEsBhVIDZd5n4f\nwqcmj8z/QUp2AgTEZAPpXmVO16Hlxize+6oH0vf9y9/LwVvxuzcaG7+xbsVnYmsv1Kb1SG3reDBi\nuyeicdbTJ8LKibexWnFJT5/Hx5bCZ0puN1LmqBRm9fWUq/vdHmiHyO8sNquYVCfmreuazczz42Tz\nqcAG6/FxNgp7XdCO00Sfab1OqNrCI+ZtRPfiK6WdoQTUGByyUvW3PP4W8KfeaqT2wdSIFeSNFY3A\n7Qc6AJ9vYKWL+yssb2LRKXYdTpYCWPmtXq+k2azGGl1rxrCZeZ6z2Ubm5Dd6nVCZaJje828veDpT\nAD/Qgq0z4+GI1887KxuBKyvNlQWwilY3bK4XU5xlnudsXsF0HdcJxS8xHjn4PQm90MGXtVpw1/nY\ndoVG8rKKIvDYC6sIzEpDmVAARh5b8E5t4qwegZzNKkIrU+kT6ncHXDz+U4fb0Eq1EzNFI3lZPREo\nd3569YBZcSTZbM4eWXHQpw0AfUJdMaC/Be/Thic/R40f5568qhEpyNsrGoHSQndF119di2ezOXtk\ndXnwI4ym9rpy7nH/rdGPsMO5a3kE8gicQhEw/0u+4RdPIcw51DwCeQROlwgc7mhP93i/JjhdvM/9\nzCOQR2BVR+Ciaxjem7mRX/MI5BHII7BKIvDzhOP/AG+s8i3bJsLDAAAAAElFTkSuQmCC\n",
       "prompt_number": 27,
       "text": [
        "\u23a7                                                                             \n",
        "\u23aa    Jxz\u22c5(J\u2093\u22c5P\u22c5Q - J_y\u22c5P\u22c5Q + Jxz\u22c5Q\u22c5R + n) - (J_z - 1)\u22c5(J_y\u22c5Q\u22c5R - J_z\u22c5Q\u22c5R - Jxz\n",
        "\u23a8P\u0307: \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
        "\u23aa                                   2                                         \n",
        "\u23a9                                Jxz  - (J\u2093 - 1)\u22c5(J_z - 1)                    \n",
        "\n",
        "                                        2        2                            \n",
        "\u22c5P\u22c5Q + l)      -J\u2093\u22c5P\u22c5R + J_z\u22c5P\u22c5R + Jxz\u22c5P  - Jxz\u22c5R  + m      Jxz\u22c5(J_y\u22c5Q\u22c5R - J_z\n",
        "\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500, Q\u0307: \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500, R\u0307: \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
        "                               J_y - 1                                        \n",
        "                                                                              \n",
        "\n",
        "                                                               \u23ab\n",
        "\u22c5Q\u22c5R - Jxz\u22c5P\u22c5Q + l) - (J\u2093 - 1)\u22c5(J\u2093\u22c5P\u22c5Q - J_y\u22c5P\u22c5Q + Jxz\u22c5Q\u22c5R + n)\u23aa\n",
        "\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u23ac\n",
        "             2                                                 \u23aa\n",
        "          Jxz  - (J\u2093 - 1)\u22c5(J_z - 1)                            \u23ad"
       ]
      }
     ],
     "prompt_number": 27
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "L_i = aircraft.linear_momentum(frame_i)\n",
      "di_L_i = L_i.diff(t, frame_w) + frame_w.ang_vel_in(frame_i).cross(L_i)\n",
      "eom_trans = di_L_i - F_T\n",
      "eom_trans"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$M \\dot{V}\\mathbf{\\hat{w}_x} + M \\left(- P \\operatorname{sin}\\left(\\alpha\\right) \\operatorname{cos}\\left(\\beta\\right) - Q \\operatorname{sin}\\left(\\alpha\\right) \\operatorname{sin}\\left(\\beta\\right) + R \\operatorname{cos}\\left(\\alpha\\right) + \\operatorname{cos}\\left(\\alpha\\right) \\dot{\\beta}\\right) V\\mathbf{\\hat{w}_y} -  M \\left(- P \\operatorname{sin}\\left(\\beta\\right) + Q \\operatorname{cos}\\left(\\beta\\right) - \\dot{\\alpha}\\right) V\\mathbf{\\hat{w}_z} + D\\mathbf{\\hat{b}_x} + Y\\mathbf{\\hat{b}_y} + L\\mathbf{\\hat{b}_z} -  M g\\mathbf{\\hat{i}_z}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAABQkAAAAmBAMAAACVEgCCAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAZqtUdu+7iTJEIs0Q\nmd0e9zMQAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAREUlEQVR4Ae1aC4xc1Xn+Z+d9Z3ZncGslcRt5\nsAhKW6IdVUkq9cHeAqqj0uBtIIE6WDs4hlbQyAt1ukoo8gSw+oqYiSpalEr4Wq1NWwnvqFEh1Qbv\nlNahLSI7iqCkDeluk6AWguJl7TX4efv953HvuXfOnbUJiPXaR9q5//nf/3f+e+651yZ66yNVfeu2\na8zy6TVWzwVUzvUXUK7vcKrFzjsc4JL7BASKtyYILkb2NRdj0auh5sPN1ZDFKsmhtJZPJ90BIJcH\nyAbZSbOVNQa4F6IPraRwMcmd02u3WueH3eTiDn49UTbQTlglaXwv0WdcUFjDuMdrHTwfqUG+abDO\nhSx9+NDWxPSLn37ZSxIOspM2CRrlU0ku+/ilXh/rImXM8atJq7tWqy8fo88l1jbn5scThAPthM3K\nGgmuQ/a0F9KXKCq5axaELiUf/iDpJhU+yE7arKSRHFfH/G9NXLoyAqm+LWFlCFcPcqs01wGHToXd\nTeeEYaGq1FJNq76QO5d59F6r2Mr8oubabYR43wZKiKiNI9eV0jSU+1URi8g5YegIcmUIg0Lipsb8\nrQBE2igJgj7gBuTa8l0jnShZ8I9FGcZskF2gdntdk+t8JjPLn8kcOVWnwpFtmFkPnanRo2yz5VOs\n77zJtB6Z0eOb7npJzfa4mk10ICA/GFAmIeTXln+HhpomW9DOfVOXa6bh02loZmBTMG4JIc7Uv9Gk\n/ojr/E2bbtPW5nWlNA3dPlURi+hm1jk/CBuh28EAJgD02Iv36Js89MSUyJGN+iEgfv+UwN3gf5Xo\noeVewnJf42+/crtHhQmXXUbG0xNLHqUnbqI91i5c5/991WpHkRYiSvudwO/oJMj7cS+LY/bPM99+\n6Bx6A6LCaJs1Uov8G4wK3N3eldMW+1MjfI8vaVbkynLcsX9OBAexcW+ThjTX8FnqBnpaWuDE1BDi\n91KxQZaIx4mmG1rTuK6U5iBVEYvo54TO+UBoFEI0CMAEgPJf01GN9AQZoNoPgXz/FMAJMl2DhX25\n8aaawr4z6wqfkZ8WO0APl6xdSOId12ZHZgthK/Nr2q2znsl1i0TZcRA84eNol69q/Lu87vsErqlR\nMRlekDz1y10y5kZYPMmArYZBapaS7yO6mmhzyJXUsw3ksS3OjSj22UCbeU6V0gthcJU+7iCANhTZ\nxKV7IzdNlj0piv9qOfiSlLHwiBCa2SQIlR8Twkj2gwCMAqQzKz8Bn4d7yrF5EYkJIyNbU0GB/Unw\nhpgfXW4VQNzaL8S7UApHThLl231dqCyLYlOwdqFZPwI/wNHFSB3ooLcON+BzEdA2qX98WLI23IXr\nT39KTEYakqd+tyClLb0IiydZDz9ylMc1ZVyzHiau+Np2GFRkHOKZbPkI3/h41WcDRRZngG6Ngogq\nffDRgQxgfGS9gKON0m7AihBZL5hKVRkL96AQJEKorEwIjUIA3gAA3QhAOrPSAnwOLSrH5iXrYeay\nka4G88iQwM2C97EIX0xUgAyfdPdWY3uhFGbQZz8BcWwvVJZ5cUa2dqFZfyTy8DD642+nXZSEBXoM\na+X7v+KfpQk8vtRQy9jDneF8BfsJRqWnZPKyE5fRuqT/pimv+H00oMQGbswkKeQuFbEgJWWtdWQp\nu5tyHvosLmiN0Mb5o4AnxOiMh70wYtCFuXF0Sy/QDQhLmnqtAx1F9KmqWNQS+SdCqMwNCM1CiAYB\nGAVIZ7a7DZ95y84uUZdGeGrahgR71CNnEjeAf3TeX8j6DaWpAuQWMZ+epNnNM091lQhHOZfJ8jIV\nmSgdu3fmJWbIoSz5VsQeGrFTGkb94FT8muLTcAZG7sYmSkIP98CeO0qtBRqpaw3SXYinTiYl6664\ngZQJnAScM7g++OHrnZ21wj98+7ormH2f+Pl4W9PlF6fq9Mx1n6fCtdd1lLz4J88faiI4RzbGNMtp\nHqZRn1LvPuGTafZ0YIn+dNuD17gwEOJHHr/+/0CL6LgGXVipUfoHYNjSIOlSG+m1tqrmn/9PeJGq\nKhbNVZmVCCELMXohhNGCBwAYA0hl5sCCaJi3neemfo/Kuza3GSlmct3KSEPAbDFC4IjQx8Nglv16\n4SzR95WGajSqNMDYW6PZ16kVbtsq+s72H7J6ya8Wlj2mxFDCSodnUTupYNSf3nVtw5mtKT7tS79J\n6e4opqlj8vNG7jjNnaD9WkEvo1Mda9Ivq2eABF2r8Bni0UkU1KARytWotIOmPQhvw9/BdnmrUyf6\nKdD764WthVcp536DaEHJM+WpF5oILlIHT43dVSZ2NuM+c8xXPtlGeNpO9ITrYCsnIXafue8W0ByR\nR9CFrSt+7XVm2NIw04SKAtSq6rxCh72/Uv5VLJK3ZSKEHBZ3qgGhyFSycS8NADAGkMosLyrJYU9w\nbqQ93kMe7RDowyOjrow0BDqOARzRRld+INvZKPttZ1zrqAAtF4wt2At7lD+lZRqXsat6zCqdQS9j\nEdRQlnsneR61kxpG/b9P0z2ar0k+0QbnKP0z8WmveIYfyPwC3V2/RFNaQS9jypuul+vZhuBPV0Mx\nal6a+u6vg1G8lfLchXj0VViB76+biT6b6Yo3EJqh9Jsjk1R48+Eucaosz8O8g6ALoI1xxOPJcp/P\nkTrYyifbCE/owtuJeF2EGK55U7gafzyCLtzSpfUeGLY0zDShogC1qmYnKdf7DeVfxaIKiAEQshC3\nmgGhyFSyBwMYA0hlJs+3lXG0QofGmv+Lf0T8Y0Yfw0BVQ6DjGMBho+vxAxmHlNOpiVpekDxVAbY0\nQeODCLop7XugxVDC7JKYcReO8c0vhxKKs0LMTiqE9adPU6VudGGPjpHrHINa+TjVhPaE+1VfHPOl\nrVrGYdzvGZp2BVMcg579EY+vYQc6rTSP7PBEF3bk3vBbeGZA9B0Gp4UASFocoN4YOc5NS5BTBmu3\niN0Af1Q+xA5/2ATp44/SnFbUZ7Yd+mQb4QldeBehDrwPQUwuEMYFEcUIuhCdWunZ0zDTZCMFqDXj\n+S4NdyaVf1fGUgdlO4QBTiaEnGkgGABgBKAws8oipzmG8HMoqXAUpU0yUhgGqgKCIEqwGAJsmmsM\nsTplj5XGFrMekzxU6TuZxi6ALiyIT8o818KcbD3uwhZCq6Es+a2nbNoFCYT1YzNreZEu3Jb25Gnv\nbArIYMy/ujD6vo4gie6dmXl5ZuYprD1le79IG5uC36oqsbhUFtTsmdkzogsbQRfy3vm5P4AYgBS5\nT8YQ43j5Xr8HmvF6BFvEuOpCTPV4mYmRBfxEffLiaZ8MpvCELrwj7MJik2g9DMUSBOmDgRt4umZP\nQ7sURvmZmUOfnplBhtaM8dzI39yV/nUs/bpmhRCh5TAhFPeLFgwAMAJQmNkIMCP6Ef72umiOE+jC\nBiOFYaAqIGCeHLpK2YWVhfcJ9vDZg0Ov/6XUCAPwy2kKjxejC0NhtiPUzS4MhWy5wbSTrvEb1o/+\nuQGNVlMip0oz/0RDoqYjj0tm69TkXnz+DobcTPZR6etV+oDkThtirEdPch2PNrZzeCKrLsSzYW4S\nR4c6xI8am9CZv6MUcBPPDsiHaljvBcyNMco0PxViPkfgS/tkG+Ep7EIW8867CX+IKIbeCwvAdKxn\nT0O71EbxvdDIuPwGQvAGwP51LPVEpiQIZSYmhJxpMAYAyEcXEyCVmTif86encC/8L0YfDH4iKyMN\ngQqkq5Rgj5zsCX7Bf7J46nKlorc7/rZKhxuiC4t9T+Rp+MfgLtwbfyKLO7fHXWjYCXUK69/o0WeN\nLkx59IIrdhe8C/yuVC75XhZ78dV+d+I15shl3EDDN5E4QII157JAD/7axSM9SUNVowtvwx1aRRey\njA/KM1DJIrsTv0DEZUNOHbwteP1vJ3vbj32scAvkMZ98qNc++e1EeAq7kMUltO5WWHJEHroLec12\nu2VrGtqlNkq7wtSaMZo5vwgx+9exRE6YxyBcv9xt8dakhglh5O1kAIBxgFRmxRPOruaeJhzjXEj1\nbWiHP2P0wTBQ1RCo+LpKBg7nqONdyZ89SRPjOAQ9IZ68MoBocBxzuJvyp9JHTmfPsraKLrpdduHu\nTlCkFPLXo0xV2mmfMg52SN1CY0S4kedrSmFE7Dgtl/Xm21I7dZzyfPzc+RWJoFzGBf7qK24RiCqT\nUlX+znblNX0r5TxshHovvAebcI1S8wXIvwuV/VX6AlqrNPmRLt2JOeTArPAKLsM9/Bgj81ovv84D\nI+YzXwt9so3wdEfwRGYx1vchroQj8tBdWFrk98IHrGlE0oSNQtua8beI3r/4Jelfx6LpuogVg9A5\nS9muEIgfE0LONBgDAIwDpDO74Tm6+2fZgfOblPf2NOlJgRQYBqoaAhVIVynBFp95WLKxQbtrVD6Z\n52epKp3BerSH6ewttGcr5U7sb4dC8XkT09Jys/BJLyhSpjYCy91daad9silGWP+cx59W52tSYWji\nB1ilj058rwuldayJUb6RHN6GRpbrYi6W8QNn3bSbvv3sVYKVrYmL/PmOv0MSxbt/9fOF0aV/HV16\nYOcnwOIHwtRl1eFduP4j/sof2uzR41M/Q3/xk+/hoiAvfnHXVBdkSRQJQo/7f/uGNjWhEPUpHibK\nJ9uwp/v9q+73P/jR5Zp8sD/yns18EBUR+aq6cHh2qUrDO1xrGpE0YaPX2pIx5e++jF6EG65Ix6KW\nzD8O4Wj7/dBSIwKhfCpKySAA4wDpzApX7nqFCmz/9PNfIOfFe9qMFM8NVDlFcxjA4fFTU6LDddrf\npHJv7FZmiABDs0ubrvw2T5967soncTlyI08ULv/hXyFmpZe+ddckv0nLIoXl07CcXVJ22qfQJqP+\n4hVTC/hq7UcVpFrs11nuCo5aRlOa65izJHrECyTh96iAhTb3xMGKOfsNtiJz9VSzn/tqyLLYEMTA\nhUcQ0ZK+1NC/ljT1WmsVfY2rqlg4cHa1SuTa6hj5RiSc6YrDApCR2Wfk9764l9AogCCmYgMOL4j+\nl7pQNAIYZtOiQ+1CXaTFUvk0HGky1xNUskKgOPqaIL+sGeE1Px7SyVRKRmKFvGtRg5yP9zx2yYv5\nm6v9kjlVtLG3WGwIYtUZQURL+lG/ljSL1aiKnsVVVSx8B9Aa0WtuqRZlhDOjkJAZoywAGZltT3di\n+mIaGgUQxNRswBF986jQNwIYZtvFsZCsQl2kRah8Go4k2aS5tqCSFEKLTUVxUAgZIYWvjucy7gyU\nng0ok7iT/kVOyxY8h32cc/tGztMsmw3OpWVXKtgjauPIdaU0DeWoqo5FhBck2yj4no3NvLCQJA3w\nBwI07zetpoFRAgRW4PAaMPnNutUfmEO1vY0kGTZRL0mW4DN/WgOWoBD6W3+MfHkaCHmaKuNF8RzG\ngUAnXL2ABeIAfw3jYVuSwss2UJyGMEiwIadRbEoFe0RtHLmulKahHFXVseQ/3BhqAXlHQMWJsJC4\nxJgPBOihhNUJjBIgsIGNI7V/yJf7kxFfkakJt6JeEvqF4tXQxk72mb7843JtBwW1u4xwE+79iA7O\nZlXFSHkxiZwG8g1WsZW5T3PtNlqcEFEbR65BGisbJagWTkYc6slfE/+3ioShM00QC3YQzl6s3VQb\nJVVzPr7sESLcgUVGNN/2Cf/jxKURIJDvBKRJ7H6wZk7XJv0uFjknd9S1iev5V5XVm37U9OBT0fma\nnL2LRWbsN/+ahPkcivrIOehcUnn7Ebjl7Xd5AXv8nws49ws59dg/UV7Ipfz4uct/mP3x/VzycJ4I\n5NzzNFjL6vu9tVzdKq6tPOAbxCpO+x1Jrfz9d8TtJacrI3DQW1nnItGI/C+ti6TmVVKm01wlibz7\nafzbu5/ChZzB/wNQ6eA+Yxy6TwAAAABJRU5ErkJggg==\n",
       "prompt_number": 28,
       "text": [
        "M\u22c5V\u0307 \u001b[94m\u001b[1mw_x\u001b[0;0m\u001b[0;0m + M\u22c5(-P\u22c5sin(\u03b1)\u22c5cos(\u03b2) - Q\u22c5sin(\u03b1)\u22c5sin(\u03b2) + R\u22c5cos(\u03b1) + cos(\u03b1)\u22c5\u03b2\u0307)\u22c5V \u001b[94m\u001b[1mw_y\u001b[0;0m\u001b[0;0m - M\u22c5(-P\u22c5sin(\u03b2) + Q\u22c5cos(\u03b2) - \u03b1\u0307)\u22c5V \u001b[94m\u001b[1mw_z\u001b[0;0m\u001b[0;0m + D \u001b[94m\u001b[1mb_x\u001b[0;0m\u001b[0;0m + Y \u001b[94m\u001b[1mb_y\u001b[0;0m\u001b[0;0m + L \u001b[94m\u001b[1mb_z\u001b[0;0m\u001b[0;0m - M\u22c5g \u001b[94m\u001b[1mi_z\u001b[0;0m\u001b[0;0m"
       ]
      }
     ],
     "prompt_number": 28
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "Lagrangian = mech.Lagrangian(frame_i, aircraft)\n",
      "Lagrangian"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\frac{J_{y} Q^{2}}{2} + \\frac{M V^{2}}{2} - g m h + \\left(\\frac{J_{x} P}{2} + \\frac{Jxz R}{2}\\right) P + \\left(\\frac{J_{z} R}{2} + \\frac{Jxz P}{2}\\right) R$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAlQAAAA0BAMAAAC+3btpAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAiavdmTLN70RmIhBU\ndrsq+tydAAAACXBIWXMAAA7EAAAOxAGVKw4bAAALAUlEQVRoBeVafWwcRxV/57vzfeyd75qGjyAs\nX1OpEoXEThOagpR6qUqLBI1PgIiQkGxEakAgEkooFaDk0gISf9lEoEhI4O3HHymI+hQQERUlR1sJ\nVWrijZEQAiU+SlBVUerLJyGJc7zZ+diZ3dnb2fj8D8wfu++9+b33ezs3Ozv79gA0bf3ZpsYaMqV2\nfjRk+z8zpJyBKZNLfhfcYQL7X8aknexVk+v7IEw6JrgQJtUImRIZrJoB/J0GmB4QIw7INYpXegQR\nXY/BfEMoSYT1ScA67Pd0RtWWslU9sablSN1+gcd98kPrZjFo6Vo48qFuC43lf79ZHr/cguz4PtQW\njNa0YKxUzbP4tIvd1zbv1oVaN75NcS6Mnt+4GRMsTClmnbKOGkWEQ92NG9/QAP0sWCdPRs+RusRg\nqYfBW3/SriboaBWNH8f5dnwKhRcJYgc5JG5/YKMiaOEywNC3dHFG6qq10obiDTT9TTWHtRIPJyKg\n23w7DAQ/C9bJk9FyDO6lKOsreD5mA+jukOLhGvYe6gDkCJwoNzfJrcfQlTROC9lzOJFXqFE9fl5V\nYdIFGK1jkq1AR1A91mAWHiGL8Qf4lJDRIgtmFMloOSo1Chscw/NABwo21ZVj4bkpnPnH2niJHbwR\nMV/4oQIwVYZshuS0UMapWr6o8z8fMC4h7biDA7s30BFUz3IDj1DGYUrrnlYiC+YhktFyzMxS2EgD\nzzghD8J3mZ90ymQwvSPzLg4mMj6JPflqvi4BTEWxwHFaGOpg0Osa//x/AsbXcZS8Md0S6AioYrEV\nEYYw+wU7ACOqyIL1+cnoOEYcD1bEuxQgcyX75Vf/6RmUQ6aMF+TO4eCkcBLY2PepV082FYiZso3D\nGC1AZQrXEZubpXMKf5QXFo+cciG/vfFp7MA58pJLAAsOOUa2tM26SASv4eQpaZdWkoXM4Sej42D3\nc8r7DYcu5bpdHJVge7Z0CUrNUTQXVqCAIwYHut0gyEAv8OSBLyMwWYW87ickq5lVLexIT8HvTt9f\nA8ieP70863GIsdAzzjnMLlaiyU1nglOUQjALhcNPRsNhsfuZ3sqVvXr2Hxevw29hP3bmL3r3nx4W\na023GYTTAizdtv4hR+d4vAapZupayYEmPIWA8jUofNMDFjTbGSnCvVzGCLQtNeGww2T5RLJQOPxk\nNBwlNtwVby4tuHIkX7ZhBdwieU5ZN6Dm25NKkw3mwWkBlpsRQRZmwYLBNunN1vBAFhy6/FskkciW\nxSWCNoxAG65yFZvJ8olkoXD4yWg4MuwXSmMiAFvlQJJsw76SU/AetxcKDakjofggx3NagAvcFDxP\nkDGsuMT8HXLABafIboFpokc18gSizYtARHSbr1GbcqRZSBxSMgGO4hZIj1HngQ6eyUNV14qzMP0b\nGPByGH9GhzC0PeLhJFrIavcJ8HSV/myTdbz9SjZxm5xlD0C8aR1iiGi5mtfBI6CSxckTfgCKLHwO\nOZkAR3mFPC7z6yDdzF8pDtd31fX0BQcedCHXJr3LmgcksZs0kjM2SgvDH16se9sqjeuSS570peYS\nlJvwsrfpxW1VagUHDtu8q3HhJhxS0miE4t3PN7wZMOJaHMDOLAuZg2yreAtwlFzrLlzYKs1dAH/8\nFWy4k+MC5zQy173dMsBEI9CZQC1f9cCUtty8I92AQf0SfQJ+WcNtT2sZnoXsBzb8qQGwhzx+sy0S\nIdfGQ1Rj04FGyLwyTLfNc+5nAg40C4VDTibXVvGnFlsAP38CzuAs3Tx8D2TVbqYNHNgBJ+DkgbvI\nT3pICzEzkpWZNI82C1/EjdyeG9oZUh5+H+KOnn5hwyxk4PAXAB7ovgdg52nPn8fxlOAB13DSaASw\nWshxHoP8K8TjZSFzKMnoOW6DTV70N+FH3nnNDrm2EprOMcVkpoi9pQ6+XzEeUbQEip7jHvi6F+O1\n0lSCWDFQ6/WxEOK4LZuyHf0kljF6ucDWFF4zkVGW8pwoNrJNuddc5hxqAef9+Y4XYqJbNw8Vixy0\nQ5AZVzLNp8ZSkppEJJUCr/GaieSc532e7al3/1TqSyIKjhtyAefxM1NekE98W431vKom1Cot7mA5\nTFoQJjQc/fPwX5ldnAwZ+cwRNRPcSTosSEl59XZ3fUkEZ0JSDhx6UcCxnKHZYDxPP6G1mhpnmhxZ\n4pMJH6O9mynjeRpG1ExwX8E5ovaGgjgpB+41RQFnoKotJQM+9VbTbhXO4jImHGHTC6aMeFOQNtQR\nBRzBkYl7WCTl2CvtX0s/qFLm4NE0KPf7/emTcPDOXyy6B880AHa/csahPeIy+jZUj9LAfs3En1V9\nGyrOEVXAoSmwY8KhSrWLF6H6lmutNLCIaF12Bmo0kBiq0aYSP6yYMn6Nuvo1E3+oBvR7Wp8rKUdU\nAcePiJJpUOa00ITr2dbdUDwHhQ5ghWuoTXv8oWLIyJMp43Yawa+ZSEOFN2XPlpSDFnC6EY1R8aAR\nKK8G+OutpD1MPPZD9qplXSFfxwbakJmCnEsDBYdKE49dH2esaCDURIFsqJabLFNpWffe+7EOEY7A\nsJwjDOAWhQM3/xVb0GiFt09Pb5+e3q3t0xrxXbg8BqUOZMbIIKWrMF9HYGp6+nOPTk/bxGfCvzSi\nBloCxq3UlddMZA6/vBMI76nJOfBhqy3gqNH5+KvWKA2LwUM2mUw4SDNO83gD3kHrAP7NMeJEOTO7\nKSNdq+SaiZi5+hcSiTchh7aAI4WjomlQisb5hLvOXBVwkB4otuYAvlqgs0hcRt/2VfTpJNdMBEff\n9lWMA7dV4QLOKofKugqjDvknAw7SZ1PkQ/UV9uIiLmOuEeJQDaY/Dt1XyTUTwcHf3dTAkpaQAxeu\ncAFHCkdF06DM8fFTtwI8AbCIX+rfBm/gd+oq7RGXsdBg0KiTISN9sVFqJoKjJJXntDQJOUIFnOL9\nb9VDgcNBtTDf76wvypK4DOV1GeCZnf+QYSiHGUFHKV5lfXfBob4u95WDsh31voX6zJ70/YCOr7g6\nGENlATcK2pafZebjttxv1WCkIRsAwoxaSvWV2AshOPirNA3cVw4a8r0A31Cz1mq9YCPNoZrWyTfm\npnwZ396a9M8isi0s6yh7PuXYCwmNtAYc+PQdDd+BocR7wV78ycdC+ICBvN76LTfmVzd8a1DSUbJt\nZhBK9WXZvAYce5pGQ2UIk5OVZfxALDV8gIlCkGQOiDrKXDsAklX8lOy3teEYlyl8sqBkCAu6eXpW\nKbuhKa2MndYHjSHKeTcKivbJltrZf478OZUhQjOERXg/ErDPVAMGnRqmXHJ0OGZjn0wFov8c6bYI\n3kswhEWEmAjY7w3oWjVM+XctjhkzHbW3/xw7VYIozRAW4T7TUDoybUWNUMKUypeGoFde3bL0n6Nc\nC1JqdUOY1heNaVvpwX19fAtTav7QI4fZJivk3SG+JeJ4GYpOfEj8+4ARLCqSepElG34WhfTtYcrB\nmt+rkebk51P/ObI2lB0NbcBkCAt4Seo+SYZPAvxF1rWyhnLG0SK5MWdzCc/953j6lvW3SwRRoiEs\nyh2/fUi/eHHLLfeNRUNZj4bybG8n8bdZhK0Bx55u90LvBLxeQ1h0pEHb7xvEAm38UIUpxT/4/VCq\ntNtX14zDp1grqYiVstW2Y42YCF6hOgYT0x3LEePfj+6XVh8k5v7DTxLxczUui1iOuAB96M/UVhuk\n0I6NcF8sIgZgwBEToR/dm1YbJOJfA3LYAVfWbkI24LiJqEldMo2kHirewgJWbPtILKInIMzxXzna\nszGZhY41AAAAAElFTkSuQmCC\n",
       "prompt_number": 29,
       "text": [
        "     2      2                                               \n",
        "J_y\u22c5Q    M\u22c5V            \u239bJ\u2093\u22c5P   Jxz\u22c5R\u239e     \u239bJ_z\u22c5R   Jxz\u22c5P\u239e  \n",
        "\u2500\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500 - g\u22c5m\u22c5h + \u239c\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u239f\u22c5P + \u239c\u2500\u2500\u2500\u2500\u2500 + \u2500\u2500\u2500\u2500\u2500\u239f\u22c5R\n",
        "  2       2             \u239d 2       2  \u23a0     \u239d  2       2  \u23a0  "
       ]
      }
     ],
     "prompt_number": 29
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "v_cm_i = point_cm.vel(frame_i).express(frame_i)\n",
      "v_cm_i.dot(frame_i.x)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left(\\left(\\operatorname{sin}\\left(\\phi\\right) \\operatorname{sin}\\left(\\psi\\right) + \\operatorname{sin}\\left(\\theta\\right) \\operatorname{cos}\\left(\\phi\\right) \\operatorname{cos}\\left(\\psi\\right)\\right) \\operatorname{sin}\\left(\\alpha\\right) + \\left(- \\operatorname{sin}\\left(\\phi\\right) \\operatorname{cos}\\left(\\psi\\right) + \\operatorname{sin}\\left(\\psi\\right) \\operatorname{sin}\\left(\\theta\\right) \\operatorname{cos}\\left(\\phi\\right)\\right) \\operatorname{sin}\\left(\\beta\\right) \\operatorname{cos}\\left(\\alpha\\right) + \\operatorname{cos}\\left(\\alpha\\right) \\operatorname{cos}\\left(\\beta\\right) \\operatorname{cos}\\left(\\phi\\right) \\operatorname{cos}\\left(\\theta\\right)\\right) V$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAABM8AAAAVBAMAAACj09SnAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMmYiu80QdonvRN2Z\nVKvu110NAAAACXBIWXMAAA7EAAAOxAGVKw4bAAALl0lEQVRoBe1abYxcZRU+M3dndnZ2ZnfBGFIR\nGFsVNCaMtX8abVhbtTEhskbWLzSMH4SaoB012iaAXCXYFbFdbaxCUEfDD0ljupGY+IFxFHUhSJxE\nEvwg6eIPjca0hVYom354nvN+3PveOffe/vIXb7j3vu97znnOc573nfuxlOil9pIC/xcFKh0iHJl2\nMcYR0asy82ZY6eBaIarOoVPQKh2i+sY8HIKZqJuLIzSIk4BHERfdbMKRYrwVV0iVDkK4xJy0sNoW\nkHSTybXSEQnypErz0H0AQNHGbh4TsbNLnh0m2yod29G1JLGXpiper492hEg2wxbOjINF/ZZc5BTN\n49IgmvE7aUNfLHIyAWymZ5JJ38t6PvDJYQrHu6FjgGIdh9WdhxPFhkcRF/ZKzIhBs+FmwOeSCr0f\nOkmJ47gakiNpQDQJNKkcR6e098kC0O31tyoVSjZDlbsaU8PHnctcxV6aKi5ar+Y7nhYihsyt5z9M\ndOmLQ6L3MgccLN5pucip1sMFO2kBHbQrYrnIyQTAXEsmfS/jWTlEz6dwvBs6BijWcXi2ByesofAo\n4AI3b8YAzYabAZ9LKvR+6CQljuNqSJ6koGgSaFI5jk5p75MBoOocPaBUKMkMVXTHFBB76lTmCnt5\nqrhovTb0ZntCRMjU/8uI7TWiSQ7CkWl7ZXwzn00vNNuAzTzrYz8WutgRzNMdOqvj2Oh6P4UTwJjk\nUYcn0S3ggjBnrncx4uYmzCg8G1teha4ulOhh9BKt2ZMM05jKRAIvVeAR8tB8ZI6J/C1hEiCkYlzB\nXoHAMVDZuQYeF5aqeL1+TLMj0cRkuJcTzPAxxYuCI9MOyvj7fP5ExoKhDdjE3fqytT9lr+EFnq8j\nek7HsUDNuRROEG9oVMEPPAq4IMyZ232MuJlw08+eiytMl+hx9RJtFk8ym8hL4KUKPEIemo+ozSUd\nSyoMEOxqYG5MgcDRL1zaNfC4sFSF61VdptmOEDFkbuQEn+bjsD34El3EJ9OaK3KN+VwbmSl645zt\nmCAexJjgR7C01CpkPB+mCt/RPI71lwuymye0x5EZe7I08FyReI+RyWDdndltNBturXwpqTBxDEt0\nuKkSx5E8SUHJEDQSOKnSebJKO58MAPWpyVvSMUkjyBLaCWd3CgR+PDB6y6xzDVzEXppKSnVE0/Eo\nZnpIV4xSq3W8S1HMTveYo3L7noUt5+jl91/yhT7PNIZ8qsw98ifbveSpz0XXrVU+8Os9/2QDgur3\n9CmeZGAAoNlVGPdsrtPkusW5ZelrVN+3d0DshiAJvqi1rWu7GTto7FrlV8VbTbzQGs/AwbcNGE3M\nfHUyuzFPcTXFFdYfXRrRrj13ix+7p0t0OLZEBSlNErVlxLISSLVKiWml4TMOQM2Xrd4wVySiroBS\nVSCWZi9NVbxeRwZ0x0CWwsh2YEAt6L/dHB8nWqF3EV3fj07yzESHT09EnRp/yC7wrpqnaZpYo9pO\nuob3BIIa7RPRaPJ5otfwGM2sguI5c+7Yv07wKybjRN+hDd1Lu7RT3DgI2VujzfQeg5O1g0Z8lLbR\nHT2Jz+FCbx7UD0YjcWFAv9GkCkygFVdIW0eVg5WHaKIvfmGJSItmN5qClCYptYViWQkglVZiSmnI\nqQDQZH3psbkCEXMUUKoKxNLsZalK1uv4sWPfM0SNbFf3zV9dfsnF8/Fgj2JstHfK+xRNj7DZqt3G\nGrVXiJqH+K7C2k3I81cCtrVWqj1alldULIJdBcVzeoFmF3j1Gae2QEfnfsofsK8AIDdk306baL/B\nydqZRjR6nF5LUyY+hwt9m+jayZ6kAKq7o6EK34orpEVqPzsdU+VZ8QtLRFo0u9HGkQKS+RLw2/yY\nBBml8cavAECtySIRcxRQqgrE0uxlqYrXi35CxB+akMzIdmQoT066j4vnY/rUZwkb7e1EL/DM1ID4\ncdWg1gpVTvD4Vzu7stF4z/RNQG8qZkZD3jQ8RrOrMO45G+OhLTi7h3w9SzQbw40bZ+bHb0yPG5ys\nnWnUWZgYN1/wyOHS4hvr75iNmBnVbTRUcdMv0N5dVmH9DMeB5mlRQpj5Eh2uLXFcq4CkIpaVAFIp\nJQZKw2dcQ5rk/WcFyCJAxBwFtKrSrpq9LFXJetFpWWG/WrR7foZLEkHBtH7X+SE22tuSjfYwb7Ta\nvFm9XTeekY027zcaHek2CN+5stHuWlz82eLihxhv3HP3iB5jA1IfYf/2Om+0ebiZ7PzyGNP1Bidr\nx05pnmT7zSZeVnw8w9Q80Ve+bDdaY3Hxhm8uLg4ZHeG+FVfYxM/r6IDolPgZXVyJkjYpUUFKk8yX\nAFJpJaaVhs94hcTf/9XlfBFzFNCqSrtq9rJUJetVP0stS1Rko9mV13NJ5uHFD7A3UXU9tdFwQ3+O\nNtM1I8KXRNSlqwd4dNqNhifeV/nVqdFLPmTMz13x3D2QmyRw/I/xEQCa7LOdqBPxXe4wD7N20Gis\n8NsX7yLE53DZHRMdYFeY0dwdDeG+FVfof9tnxM/o4kp0uPaOpiClSeZLkFdiSmnIoADwb41m1kyF\nWZHw/pGjgFZV2lWzl6UqWa/6utyL3Grxk+TkUNbAfgz8mejvqY2G1+j38bvTtea1rh3TTCe10RD0\nW351+itf8YaLZlZB8ZwaNLAB8GrI72A0up9/1j8CIE8yUK1f7U4MDU7WDhqt+WqvPUheL5UMs+x2\ngIwLX/xGQ7hvxRXyOxq1p4bUXBc/85niSsx8DChIaZIKQSsBpNJKTCkNORUAYuG2dvNFzFNAqSpw\nVexlqcrW6z6+BSWrxa/6p3qyBvv5zMdferQDz03+7wWeaazhFhNjA7T4aB+iiS7fztwdDUEPUr86\n4Osf+ECzG23cs9WHRIITfYMa3Q1z9H4B5EkGaq40sJ8FJ2sHjebyJOHmCx45XFprVP1jpSdm9vQb\nDeG+FVdIWzt0Z+XrVIvFT5j5EpEWzZSoaZUmqYhlJYBUWokppSGDAsC/yMrPORpMsggQMU8BparA\nVbGXpSpbr9XooCVqZGut85DbYXP88MqrBpef/8fl55958sU1c4uu/Gb7F9lYG/DwS5+5u3L83EeO\nn7vsuu/yHILaT2zn7weiD+LEzayC4llf/T3stQGfPrV6J0WP7h/AjYcC9MpXXzXHXcHJ2OWZ9fST\nV2Kj1gb5XJY2dlr7bAp2dY9O98jjKaLiCqn+n71desvSv8WP3dMl1gY8wc1uNAWJUiTzJcgrMaU0\nZFAAmhfvW+oxg9qATxmRhGqOAlpVaVfFXpqqZL1uua3HHGsDtxzVNZI2zYuII9Me4nHUweRWnDJN\nAiojzEbL1mZXQfOUOQ3HZJb/n+hxAgDQoMt8vIaRIu/MbqORhAeAyaC4QsPMlOgl0Es0WRKSSQrT\n8wTzS0yU1nwYgD+ppbkK7dBcfIJksbwCgWOw0nlQ5anK1gs5gR5mqA75gcpHpr2Bx1XZf3ynGGsS\nMNnDfKNvrT+w1/CSQGs4JvMmhHicIB406PMyhXgNI8ngzc2ORBBJuO1nL8UVGmamRI+rl2iyJCSz\niTzB/BITpTUfBnCrrwmQXj5n9wpkyHguupaAKk9Vtl7IqazWDp7GEbYJ3mSSsr4QGswIAYbRTZo5\nNeegdRzJHMNdxwENEjvidQyXQTFLONCVVlKhMDMl6mkTyJBkMu96jmB+iYnSqs8O/HkHLYeJS5Bn\nN8FyLnO9kFQxkFSi/CYvWYRIliz+qZv/53DGjc/RPP8pFSMXin7SECBmZY8mXug5aB1HzH34eQUw\n8A00SOyI1zFcBsUs4R4s7JRUKMxMiXraBC0kmcy7niOYX2KitOqzBX/bQsth4hLk2U2wnMtcLyRV\n2XpZolmyeD8w7wgpPsR/QzNtUzDrBj7A3PTdtHJ1njqOz5yH42jwn1pIjvEULoNm9uHjYSUVemZ5\naVOIPotepCVYViKCdZ+iCkGjzJ6iWuZ6oXadqN8246v1Py9bGTXbhgtuAAAAAElFTkSuQmCC\n",
       "prompt_number": 30,
       "text": [
        "((sin(\u03c6)\u22c5sin(\u03c8) + sin(\u03b8)\u22c5cos(\u03c6)\u22c5cos(\u03c8))\u22c5sin(\u03b1) + (-sin(\u03c6)\u22c5cos(\u03c8) + sin(\u03c8)\u22c5sin(\n",
        "\u03b8)\u22c5cos(\u03c6))\u22c5sin(\u03b2)\u22c5cos(\u03b1) + cos(\u03b1)\u22c5cos(\u03b2)\u22c5cos(\u03c6)\u22c5cos(\u03b8))\u22c5V"
       ]
      }
     ],
     "prompt_number": 30
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x_vect = sympy.Matrix([P, Q, R, phi, theta, psi, V, alpha, beta, x, y, z])\n",
      "u_vect = sympy.Matrix([L, D, Y, l, m, n])\n",
      "\n",
      "State = collections.namedtuple('State', ['phi', 'theta', 'psi', 'P', 'Q', 'R', 'V', 'alpha', 'beta', 'x', 'y', 'z'])\n",
      "Input = collections.namedtuple('Input', ['L', 'D', 'Y', 'l', 'm', 'n'])\n",
      "\n",
      "remove_t_sub = {var[0](t):var[0] for var in x_vect.tolist() + u_vect.tolist()}\n",
      "\n",
      "rhs = sympy.Matrix([\n",
      "    rot_sol[P(t).diff(t)],\n",
      "    rot_sol[Q(t).diff(t)],\n",
      "    rot_sol[R(t).diff(t)],\n",
      "    euler_sol[phi(t).diff(t)],\n",
      "    euler_sol[theta(t).diff(t)],\n",
      "    euler_sol[psi(t).diff(t)],\n",
      "    sympy.solve(eom_trans.dot(frame_w.x), V(t).diff(t))[0],\n",
      "    sympy.solve(eom_trans.dot(frame_w.z), alpha(t).diff(t))[0],\n",
      "    sympy.solve(eom_trans.dot(frame_w.y), beta(t).diff(t))[0],\n",
      "    v_cm_i.dot(frame_i.x),\n",
      "    v_cm_i.dot(frame_i.y),\n",
      "    v_cm_i.dot(frame_i.z),\n",
      "    ]).subs(remove_t_sub)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 31
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print rhs"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Matrix([[(Jxz*(J_x*P*Q - J_y*P*Q + Jxz*Q*R + n) - (J_z - 1)*(J_y*Q*R - J_z*Q*R - Jxz*P*Q + l))/(Jxz**2 - (J_x - 1)*(J_z - 1))], [(-J_x*P*R + J_z*P*R + Jxz*P**2 - Jxz*R**2 + m)/(J_y - 1)], [(Jxz*(J_y*Q*R - J_z*Q*R - Jxz*P*Q + l) - (J_x - 1)*(J_x*P*Q - J_y*P*Q + Jxz*Q*R + n))/(Jxz**2 - (J_x - 1)*(J_z - 1))], [(Q*sin(psi) + R*cos(psi))/cos(theta)], [Q*cos(psi) - R*sin(psi)], [P + Q*sin(psi)*tan(theta) + R*cos(psi)*tan(theta)], [(-D*cos(alpha)*cos(beta) - L*sin(alpha) + M*g*(sin(alpha)*cos(psi)*cos(theta) + sin(beta)*sin(psi)*cos(alpha)*cos(theta) - sin(theta)*cos(alpha)*cos(beta)) - Y*sin(beta)*cos(alpha))/M], [(D*sin(alpha)*cos(beta) - L*cos(alpha) + M*V*(-P*sin(beta) + Q*cos(beta)) + M*g*(-sin(alpha)*sin(beta)*sin(psi)*cos(theta) + sin(alpha)*sin(theta)*cos(beta) + cos(alpha)*cos(psi)*cos(theta)) + Y*sin(alpha)*sin(beta))/(M*V)], [(D*sin(beta) - M*R*V*cos(alpha) + M*V*(P*cos(beta) + Q*sin(beta))*sin(alpha) + M*g*(sin(beta)*sin(theta) + sin(psi)*cos(beta)*cos(theta)) - Y*cos(beta))/(M*V*cos(alpha))], [V*((sin(phi)*sin(psi) + sin(theta)*cos(phi)*cos(psi))*sin(alpha) + (-sin(phi)*cos(psi) + sin(psi)*sin(theta)*cos(phi))*sin(beta)*cos(alpha) + cos(alpha)*cos(beta)*cos(phi)*cos(theta))], [V*((sin(phi)*sin(psi)*sin(theta) + cos(phi)*cos(psi))*sin(beta)*cos(alpha) + (sin(phi)*sin(theta)*cos(psi) - sin(psi)*cos(phi))*sin(alpha) + sin(phi)*cos(alpha)*cos(beta)*cos(theta))], [V*(sin(alpha)*cos(psi)*cos(theta) + sin(beta)*sin(psi)*cos(alpha)*cos(theta) - sin(theta)*cos(alpha)*cos(beta))]])\n"
       ]
      }
     ],
     "prompt_number": 45
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "const = {J_x:0.1, J_y:0.1, J_z:0.1, Jxz:0.1, M:1, g:9.8}\n",
      "x0 = [0,0,0,0,0,0,10,0,0,0,0,0]\n",
      "u0 = [0,0,0,0,0,0]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 56
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "sympy_utils.save_sympy_expr(\n",
      "    {'rhs': rhs, 't': t, 'x_vect': x_vect, 'u_vect': u_vect,\n",
      "     'constants': const},\n",
      "    os.path.join('save', 'wind_kinematics.sympy'))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 92
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "wind_kinematics = sympy_utils.load_sympy_expr(\n",
      "    os.path.join('save', 'wind_kinematics.sympy'))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 93
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "f_sim, f_jac = sympy_utils.rhs_to_scipy_ode(**wind_kinematics)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 94
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%timeit\n",
      "f_sim(1, x0, u0)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "10000 loops, best of 3: 21.4 \u00b5s per loop\n"
       ]
      }
     ],
     "prompt_number": 95
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "ode = scipy.integrate.ode(f_sim, f_jac)\n",
      "ode.set_initial_value(x0)\n",
      "\n",
      "dt = 0.01\n",
      "Data = collections.namedtuple('Data', ['t', 'x'])\n",
      "\n",
      "data = Data([], [])\n",
      "u = u0\n",
      "\n",
      "while ode.t < 10:\n",
      "    data.x.append(ode.y)\n",
      "    data.t.append(ode.t)\n",
      "    state = State(*ode.y)\n",
      "    state.P\n",
      "    u = [9.7,1,1,1,1,10]\n",
      "    ode.set_f_params(u)\n",
      "    ode.set_jac_params(u)\n",
      "    ode.integrate(ode.t + dt)\n",
      "\n",
      "hist = State(*np.array(data.x).T)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 96
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.plot(data.t, np.array([hist.x, hist.y, hist.z]).T);"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXkAAAEACAYAAABWLgY0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl4VdW9//H3N2GGMEuAgAQBKyAqKIIoGKs4tqK21dba\nqrXWans7t0L7a6W997YOtd5Otre2erVWhSpqbZ1ASVWQQUEBAQElQEASxjBlzvf3xzoZiBEy75Nz\nPq/n2c/eZ59zcpZ55HNWvnvttczdERGRxJQSdQNERKTlKORFRBKYQl5EJIEp5EVEEphCXkQkgSnk\nRUQSWL1C3szuN7M8M1tZx3PfNbMKM+td49wMM1tvZmvN7PzmbLCIiNRffXvyDwAX1j5pZoOBqcCm\nGudGAVcBo2LvudfM9BeDiEgE6hW+7v4qsKeOp34F/KDWuWnAo+5e6u45wAbg9KY0UkREGqfRPWwz\nmwbkuvuKWk8NBHJrPM4FMhr7OSIi0njtGvMmM+sC/JBQqqk6fYS3aO4EEZEINCrkgWFAJvC2mQEM\nAt40swnAVmBwjdcOip07jJkp+EVEGsHdj9SpPkyjyjXuvtLd0919qLsPJZRkxrl7HvAP4LNm1sHM\nhgIjgCUf8XO0uXPbbbdF3oZ42fS70O9Cv4sjbw1V3yGUjwILgePNbIuZXV87r2sE92pgNrAaeA64\nxRvTMhERabJ6lWvc/XNHef64Wo9/Dvy8Ce0SEZFmoPHrcSArKyvqJsQN/S6q6XdRTb+LxrOoKilm\npiqOiEgDmRne0hdeRUSkbVDIi4gkMIW8iEgCU8iLiCQwhbyISAJTyIuIJDCFvIhIAlPIi4gkMIW8\niEhbUVjY4Lco5EVE2oJ9++Ciixr8NoW8iEi827kTzj0XTjihwW9VyIuIxLOtW+Hss+G88+APf2jw\n2xXyIiLx6r33YPJk+OIX4Re/AKv3vGRVFPIiIvFo5UqYMgVuvTVsjVTflaHuN7M8M1tZ49xdZrbG\nzN42szlm1qPGczPMbL2ZrTWz8xvdOhGRZLRoUSjP/PKXcNNNTfpR9e3JPwBcWOvci8Bodz8ZWAfM\nADCzUcBVwKjYe+41M/3FICJSH/PmwSc/CfffD5874qJ89VKv8HX3V4E9tc7NdfeK2MPFwKDY8TTg\nUXcvdfccYANwepNbKiKS6J58Eq6+Gp54Ai65pFl+ZHP1sL8EPBs7Hgjk1nguF8hops8REUlMDz4I\nt9wCzz0XavHNpF4LeR+Jmf0IKHH3R47wsjrX+Zs5c2bVcVZWltZxFJHk9JvfhPr7/PkfGgufnZ1N\ndnZ2o390vdd4NbNM4Bl3H1Pj3HXAjcC57l4UOzcdwN1vjz1+HrjN3RfX+nla41VEkps7zJwJjzwS\navFDhhz1La22xquZXQh8H5hWGfAx/wA+a2YdzGwoMAJY0tjPERFJSOXlcPPN8M9/wmuv1SvgG6Ne\n5RozexQ4G+hrZluA2wijaToAcy0M0H/d3W9x99VmNhtYDZQBt6jLLiJSQ1ERfP7zUFAQSjTdu7fY\nR9W7XNPsH6xyjYgko4ICmDYN0tPhoYegY8cGvb3VyjUiItJAH3wQ5qEZMwYefbTBAd8YCnkRkdaw\nfj2ceSZ85jNhNE1K68Rvk4dQiojIUbz5JnziE/Czn8GNN7bqRyvkRURa0rx5YXqCP/0JLr+81T9e\n5RoRkZYya1b1NAURBDyoJy8i0jJ+8xu4887Qkz/ppMiaoZAXEWlOFRXwve+FOWheew0yMyNtjkJe\nRKS5FBWFVZy2b4cFC6B376hbpJq8iEiz2L0bpk4NS/S9+GJcBDwo5EVEmi4nJ4yBnzAh3OTUqVPU\nLaqikBcRaYo33wwBf8stYbrgVrrJqb5UkxcRaaxnn4Vrr41sDHx9xNdXjohIW3HfffClL8E//hG3\nAQ/qyYuINIw7/OQnofb+yitw/PFRt+iIFPIiIvVVUhLmnlmzBhYuhH79om7RUSnkRUTqY88e+NSn\noFu3sNBH165Rt6he6lWTN7P7zSzPzFbWONfbzOaa2Toze9HMetZ4boaZrTeztWZ2fks0XESk1WzY\nAGecASefDE8+2WYCHup/4fUB4MJa56YDc939eOCl2GPMbBRwFTAq9p57zUwXeEWkbXr1VTjrLPjW\nt+CeeyA1NeoWNUi9wtfdXwX21Dp9KfBg7PhB4LLY8TTgUXcvdfccYANwetObKiLSyh5+OJRoHnoI\nvvrVqFvTKE2pyae7e17sOA9Ijx0PBBbVeF0ukNGEzxERaV0VFXDbbSHk58+H0aOjblGjNcuFV3d3\nMzvSqtx1Pjdz5syq46ysLLKyspqjOSIijVdYCNdfD5s3w+LFkY+gyc7OJjs7u9HvN/cjZXONF5pl\nAs+4+5jY47VAlrtvN7MBwHx3P8HMpgO4++2x1z0P3Obui2v9PK/vZ4uItIq8PLjssjA98AMPxNUc\nNJXMDHe3+r6+KRdE/wFcGzu+FniqxvnPmlkHMxsKjACWNOFzRERa3jvvwMSJcP758MgjcRnwjVGv\nco2ZPQqcDfQ1sy3AT4DbgdlmdgOQA1wJ4O6rzWw2sBooA25Rl11E4toLL8AXvgC/+hVcc03UrWlW\n9S7XNPsHq1wjIlFzh9/9Dn7+c/j738NQyTjX0HKN7ngVkeRUUgJf/zq8/nqYomDo0Khb1CIU8iKS\nfHbuDOPfe/YMAZ+WFnWLWozuRBWR5LJyJZx+eljo48knEzrgQT15EUkmTz8NX/4y/PrXcPXVUbem\nVSjkRSTxucPtt8Pvfw//+lfoyScJhbyIJLbCQrjhBli/PtzBmpFcs6yoJi8iiWvrVpgyJRy/8krS\nBTwo5EUkUS1cCBMmhPVX//Y36Nw56hZFQuUaEUks7vC//xvWYX3gAbjkkqhbFCmFvIgkjqIi+NrX\nQu19wQIYMSLqFkVO5RoRSQxbtoT6+759sGiRAj5GIS8ibd+//x3q75/+NMyeHRbbFkDlGhFpy9zh\nt78NE4z99a8wdWrULYo7CnkRaZsOHYKbboJVq8IkYwk6wVhTqVwjIm1PTk6YFtg9XGBVwH8khbyI\ntC3z5oUVnK69NpRounSJukVxrckhb2YzzOwdM1tpZo+YWUcz621mc81snZm9aGY9m6OxIpLEKirC\n/DNf/CI89hh885tg9V47I2k1aWWo2OLeLwMj3b3YzGYBzwKjgZ3ufqeZ3Qr0cvfptd6rlaFEpH72\n7Anhvns3zJoFgwZF3aLItOZC3gD7gFKgi5m1A7oA24BLgQdjr3kQuKyJnyMiyWrZMjj1VBg+HLKz\nkzrgG6NJIe/uu4G7gc2EcN/r7nOBdHfPi70sD0hvUitFJPm4w5/+BBdcAHfcAffcA+3bR92qNqdJ\nQyjNbBjwLSATKAD+bmaHLXXu7m5mddZlZs6cWXWclZVFVlZWU5ojIoni0CG4+WZ480147TX42Mei\nblFksrOzyc7ObvT7m1qTvwqY6u5fjj3+AjAR+DhwjrtvN7MBwHx3P6HWe1WTF5EPW7cu3Ll60klh\norGuXaNuUVxp7Zr8WmCimXU2MwPOA1YDzwDXxl5zLfBUEz9HRJLBE0+EtVdvuSUMj1TAN1mTyjXu\n/raZPQS8AVQAy4A/AWnAbDO7AcgBrmxiO0UkkRUXww9+ENZgffZZGD8+6hYljCaVa5r0wSrXiAjA\nhg1w1VVw7LHwl79A795RtyiutXa5RkSk8R57DM44A667DubMUcC3AE1QJiKt79ChcMdqdja88AKM\nGxd1ixKWevIi0rpWr4bTT4eDB8MQSQV8i1LIi0jrcIf77w+rN33722Fx7e7do25VwlO5RkRa3v79\n4eam5cvDKk6jR0fdoqShnryItKzly8PcM507w9KlCvhWppAXkZZRuTTf+efDT38K992nud8joHKN\niDS//Hy4/vqwX7gQRoyIukVJSz15EWlezz0Hp5wCJ5+sgI8D6smLSPMoKoJbb4Unn4RHHgHNKhsX\n1JMXkaZ7550w9n3bNnjrLQV8HFHIi0jjucPvfx9C/VvfgtmzNTVBnFG5RkQaJz8fbrgBtm9X7T2O\nqScvIg335JPhwuqJJ8KCBQr4OKaevIjU3969YWKxBQvCAh+TJkXdIjmKJvfkzaynmT1uZmvMbLWZ\nTTCz3mY218zWmdmLZtazORorIhGaNy8syde1a7i4qoBvE5q8aIiZPQj8293vN7N2QFfgR8BOd7/T\nzG4Fern79Frv06IhIm3BoUNhaORTT8Gf/wwXXBB1i5Jaqy4aYmY9gMnufj+Au5e5ewFwKfBg7GUP\nApc15XNEJCKLF8PYsbBnD6xYoYBvg5parhkK7DCzB8xsmZndZ2ZdgXR3z4u9Jg9Ib+LniEhrKimB\n//f/YNo0+O//hocfhl69om6VNEJTQ74dMA64193HAQeBw8oysZqM6jIibcWqVTBhArz9dqi9f/rT\nUbdImqCpo2tygVx3Xxp7/DgwA9huZv3dfbuZDQDy63rzzJkzq46zsrLI0l1yItEpL4e774a77oI7\n7ggTjFm9S7/SQrKzs8nOzm70+5vjwusrwJfdfZ2ZzQQq5xLd5e53mNl0oKcuvIrEsTVrQqh37gwP\nPACZmVG3SD5CQy+8NkfInwz8GegAvAdcD6QCs4FjgRzgSnffW+t9CnmRqJWVwS9/GXrwP/sZ3HQT\npOgeyXjW6iHfWAp5kYitXBl67716haGRQ4ZE3SKph1YdQikibVBpKfznf8LHPw5f/Sq8+KICPoFp\nWgORZPLWW6H3PmAALFsGgwdH3SJpYerJiySDoiL48Y/Deqvf/Cb8618K+CShnrxIosvODhdUR48O\nPfmBA6NukbQihbxIotq1C77//TCx2G9/G+5elaSjco1IonEP0xCMHg1paWFpPgV80lJPXiSRvPce\n3HxzWLXpmWdg/PioWyQRU8iLJILS0nBD0y9/CdOnh/VW2+mfd0soKYGDB8O17OLixu9LS8NMEg3d\nGkr/F4i0dYsWwVe+AhkZsHQpDB0adYviSkUF7N8PBQVhYavKfeXxvn1w4EAI7gMHDj+u65x7WDel\nc2fo2BE6dTp8X9e52vtevaB9e0hNbfj27LMN++/XHa8ibdWePWE64Dlz4J574KqrEnpCsZKScC15\n1y7YubPufc3wrhniXbtCjx7Qs+eH92lpYevWLbyuW7cPH9d83KFDtL+Hht7xqp68SFvjDn/9a1it\nadq0cGG1d++oW9VgpaXh0kFeHmzfHvaV244dHw7wwkLo0ydsffsevs/ICCsT9uxZvVWGePfuoQec\nrBTyIm3JqlXwta+FusHTT8Ppp0fdog85cAC2boXc3LDftu3DIb59e+hh9+0L6enQv3/Yp6fDoEEw\nblx1gFeGeffuCf2HSotRuUakLdi/H376U3jwwbC/6aZW7566h151ZXhX7mse5+aGHnpGRtgGDQr3\nXtUM8cpQ79NHE142hso1IonEHR5/HL7znTCh2KpVISVbQGEhbN4MOTlh27gRNm06vEfetWsI7soA\nz8iASZOqjzMywkVF9bjjh3ryIvFq/Xr4+tdDwt57L0yZ0qQfV1R0eIhXBnnl8Z49YTqbzMwwQCcz\nM0xOOWhQdY+8c+em/kdJU2k+eZG2bv9++PnP4b77YMYM+MY3wni7eti1K3w3bNhQva8M8127qkO8\n5lYZ6AMGqHzSFkRSrjGzVOANwnqvnzSz3sAsYAgfsTKUiNRSURGmI5gxA847D1asqHMysb17Q4DX\ntZWXw4gR1dv558Nxx4UQHzgwuUeZJKtm6cmb2XeAU4E0d7/UzO4Edrr7nWZ2K9BLa7yKHMGiRWEK\nYIBf/5r9oyd+ZJAXFh4e5DW3Y45RPTzRRbHG6yDg/4D/Br4T68mvBc529zwz6w9ku/sJtd6nkJek\nt2vlNsq+N51ui19izvhf8FDFNbyzJoWCAhg2rO4g799fQZ7MoijX3AN8H+he41y6u+fFjvOAlhkO\nINIGuMMHH8Dq1WFbswY2rCoia9mvuLnwbl4YeBNLpq1l2ClpfHckjBoVLnSqPi7NoUkhb2afAPLd\nfbmZZdX1Gnd3M1OXXRKeexhmuGLF4YG+enWYr2TkSBg9soLLDz3K5PU/wqaMpePvlvLFYcfxxagb\nLwmrqT35ScClZnYx0AnobmZ/BfLMrL+7bzezAUB+XW+eOXNm1XFWVhZZWVlNbI5I6zh0KMwmsGLF\n4Vtqari9fvTocDPqddeFcO/bF5g/PyzikZICsx9q8pBISQ7Z2dlkZ2c3+v3NNoTSzM4Gvheryd8J\n7HL3O8xsOtBTF16lLXIPNwMtXw5vv10d5ps3wwknhECvudV5n9Lq1fCDH4T9L34Bn/mMajHSaJGN\nk4+F/Hdjo2t6A7OBY/mIIZQKeYk3FRVhXPny5bBsWfU+NRXGjoVTTglBfvLJcPzx9Ri6/sEHcNtt\n8NRTYVjkLbeEuo1IE+hmKJF6KC0NHeuagf7222E+lbFjwzZuXNgPGNDA0SwHDoTFO377W/jSl+CH\nPwz3+os0A81dI1JLYWEI8MowX748BHxmZnWYX3ZZ6Kk3acbesjK4/36YOTPMM/Pmm+FDRCKknrwk\nFPdQclm0CBYvDvs1a0L9vLJnPm4cjBkTJttqtg996in40Y/CIPa77oJTT22mHy5yOJVrJKns2QNL\nllQH+pIlIbwnToQJE8J+7NgWnFjr5ZdDvb24OMw3c9FFulNJWpRCXhJWWVmYaXfRouqeem4unHZa\ndaBPmBBq6C3ujTdCuG/cCP/1X3DllRoxI61CIS8JY+vW6h764sWhpj54cAjzykAfPRrateaVpbVr\nw7qqixbBj38cLqzWc4ZIkeagkJc26dChEOI1a+mFhYeXXcaPD2t2RmLz5rAi0zPPwPe+F+Z579Il\nosZIMlPIS5uweTO89hosWBACfe3a0CuvWXYZNiwOyttbtoQbmGbNgq9+NdyxGtk3jYiGUEocKi8P\nUwC89lr1VlwMZ50FZ54Jn/98C18cbYzNm0O4z54NX/kKvPtubG4CkbZFIS/Nrrg4lFwqA33hwnC7\n/1lnwdSpoeoxfHgc9NLronCXBKOQlyarDPXs7LAtWRIm5ZoyBW68Ef7v/6Bfv4gbeTSbNsHttyvc\nJeEo5KXB6gr1UaMgKytckzzrLOje/Sg/JF688w7ceSf8858Kd0lIuvAqR3WkUM/KamOhXun110PP\nffHisFD2zTdrfhlpEzS6RposIUMdwvQDzz8fwn3LljBS5rrr4uyKr8iRKeSlwdxh5UqYOxdefDFc\nKB05MgFCvVJZWai133FHeDx9epjTvVXvohJpHgp5qZft22HevBDqc+eG+V6mToXzz4dzzkmQoeCF\nhfDAA2Ha38GDQ7hfeGGcDusRqR+FvNSpsDAMZ3zxxbBt3hxmw506NWzDhkXdwma0Ywfcey/84Q9h\nDb7p02HSpKhbJdIsWvVmKDMbDDwE9AMc+JO7/ya2MtQsYAgfsTKUtKzKEkxlqL/+eljRaOpU+OMf\nwxQBCVetePdd+NWvQmnmM58Ja6qOHBl1q0Qi1aSevJn1B/q7+1tm1g14E7gMuB7Y6e53mtmtQC+t\n8drytm+vrqvPmwfduoXyy9SpoQTTo0fULWwB7vDKK6Eks2RJGCVzyy1tYGC+SONEWq4xs6eA38W2\ns909L/ZFkO3uJ9R6rUK+iQoL4dVXq3vrubmHl2COOy7qFrag0lJ4/HG4+27Yvx+++134whc0UkYS\nXpQLeWcC/wZOBDa7e6/YeQN2Vz6u8XqFfAO5w4oV1aG+aFFYsq6yt37aaQlYgqlt9274y1/gd7+D\noUNDuF9yieZyl6QRyQRlsVLNE8A33X2/1Ri94O5uZnWm+cyZM6uOs7KyyMrKao7mJJSdO+GFF8Lw\n7rlzw1DG888PM90+8UQbH9rYEG+/HRbGfuIJuPTSsD/ttKhbJdLisrOzyc7ObvT7m9yTN7P2wD+B\n59z9f2Ln1gJZ7r7dzAYA81WuqZ/y8rDo0HPPhW3t2lCCufDCEO5Dh0bdwlZUWgpPPhnCfePGUG+/\n8UbV2yWptWq5JlaKeRDY5e7frnH+zti5O8xsOtBTF14/2o4dobf+3HOhDJOeHpYKveiicCNShw5R\nt7CVbdsWSjJ//COMGBH+bLnssiSoRYkcXWuH/FnAK8AKwhBKgBnAEmA2cCwfMYQymUO+oiL01p99\n9vDe+kUXhR77scdG3cIIlJWFX8if/xyuJl95JXzta3DSSVG3TCSu6GaoOHXgQKipP/MM/Otf0KdP\nuF6YtL31Su+9B/ffH+YjzsyEL385jHHv1i3qlonEJYV8HMnJCTPY/vOfYT6YiRPhE58IW0IPbzya\nQ4fgqadCSWbFijD08YYbwvp/InJECvkIlZSENUuffTb01nfuDD31T34yXDRNmpEwdSkvh5dfhocf\nhqefDt94118fau0dO0bdOpE2QyHfynbtCj31Z54Jd5kefzxcfHEoxZx6apIP33aH5ctDsD/6KAwa\nBNdcA1ddBf37R906kTZJC3m3gvfeC53Rp5+Gt96Cc88NQ7d///swMibp5eTAI4+EcC8qCsGenQ0f\n+1jULRNJOurJ14N7GA3z9NOhlLxzZyjBTJsWAl530hMW4fj738PkYBs2hNEx11wDZ5yhqX1FmpHK\nNc2kvDxMzTtnTrgfp3NnuPzyEOwTJiR5GabSBx9UB/vq1aG+ftVVYTxo+/ZRt04kISnkm6C4OFwb\nnDMn9NoHDYIrrgjbyJHqkAKQnx+mFJg1K0w18MlPhmCfOjWJx4GKtB6FfAMdPBjmhZkzJ4yKGT06\nhPrllyfZFAJHsnNn+HNm1qxQt7r44hDsF1wAnTpF3TqRpKKQr4c9e8KImDlz4KWXwmi+K64IpZgB\nAyJpUvzZvj38OTNnTpju8oILQp394ouhS5eoWyeStBTyH6FmZr3+eigbX3FFuDGpd+9Wa0Z8y8kJ\nPfY5c2DVqjDI//LLw153oIrEBYV8DXVl1hVXhPlhlFmEYUNr1oRf0Jw5YdWRadNCsJ97rm5SEolD\nSR/y778fSsePPx4Wq542LQS7MiumqAj+/e9wS+6zz4bbdC+/PPySzjxTMz2KxLmkDPnNm8Movlmz\nYNMm+NSnQvl48mRlFhDGsD/3XAj27GwYMybcknvJJeFYw4ZE2oykCflt28IQ7VmzYN260Bm98sqw\nYHXSB/uOHTB/friq/PLLUFAQhjheckm4gNqnT9QtFJFGSuiQ37kzBPtjj8HKlWEqgauugvPOS+J7\nb9xDT33RojDV5csvhz9tpkwJV5fPPTeMC9XdWyIJIeFCvqgoTP718MOhlHzRRfDZzybxEO0DB8Kk\nX4sWhe3118PtuWecEcaCnnNOmBkt6f+cEUlMcRPyZnYh8D9AKvBnd7+j1vMfGfIVFWFxoIcfDjdX\njhsXpkG54ookmq63vDxcYFi1KtxZWrlt3Qonnlgd6hMnhsU2VFcXSQpxEfJmlgq8C5wHbAWWAp9z\n9zU1XvOhkF+zJgT73/4GaWlhLYmrrw7TCySkwsIQ2jk5sH794VtOTliwevRoOPnk6u3449VLF0li\n8TLV8OnABnfPiTXqMWAasKb2C/PyQo39r38NF1OvvjrctHTyyS3UspZUUQF794ZJ5uva8vJCqOfm\nhv3BgzBwIAwZEhasHj48DAkaMQKGDdP0liLSZC3Vk/80cIG73xh7fA0wwd3/o8Zr/NuTf0zuFjjh\neOOUk+C4TKOdGSlmGEaKh2+tFMLjVA/Ptbd2tLPUD22pKamkYOFi5NG2ioqweHRpadgf7bi4OCxb\nd/Bg2Nd1XFQEPXqE0Ss1t969wz49HTIywp8mGRnQt6/KLCLSIPHSk6/XN8fuNXeTlgpbV8AxG9vT\nL609FYCbhz1hD05F7LiCCsqpoNyr92WVj70ifEGkpJBiqVhKCpaSSkpKCmappKSmkmIppKSkQmpq\nKHu0b09qh46ktO9Iuw6dwtaxM+07dqZd19hxh66065xO+7QetE/rQYe0nnTo3ovO3fvQuUcfrGtX\n6No1XAnWKBYRaUbZ2dlkZ2c3+v0t1ZOfCMx09wtjj2cAFTUvvrbEHa/uTllFGcXlxZSUl1RtxWWH\nPy4pL6l6TVFZEQdLDnKg5EDVdrC0+vGh0kN1bgdLD3Kw5CAl5SV079idHp16hH3HHqR1TCOtQxrd\nOnQjrUNa1ePuHbvTs1NPenbqSa/OvaqOu3fsTsfUjph69SJyFPFy4bUd4cLrucA2YAn1uPDaFpWU\nl7C/eD8FxQUUFBVQUFzAgZID7C/eH/Yl+6seFxSH5/cU7mFv0V72Fu1lT9Ee9hfvp8Ir6N6x+2Fb\n1RdCp151fjlUnu/VuRdpHdL0JSGSBOIi5GMNuYjqIZR/cfdf1Ho+IUK+uRSXFbO/ZD/7ivexr3hf\n1RdG5RfCnqJa+xpfFLsLd1NcXkzvzr3p3bk3fTr3CfsufejdKbavcb7quc696dq+q74cRNqQuAn5\no36wQr5ZFZcVs6doD7sO7WJ34W52FYb97sLdHzpX87myirLDvgD6de1Hv679SO+aXn3crfq4R8ce\n+lIQiZBCXhqksLSw6sthV+EudhzcQf7BfPIO5pF/MP9Dx0VlRXV/EdTxpXBMl2Non5qs802ItAyF\nvLSoorKiqsDPP5hP3oG6vwzyD+az89BO0jqmMaDbAAZ1H0RGWgYZ3TPISMsIj2PHfbv01V8HIvWk\nkJe4UeEV7C7czbb929i6bytb928ld1/u4cf7t3Kw5CAD0mp8EcS+DGp+MQxMG0iHVC0ULqKQlzan\nsLSQbfu3VYV+7S+Brfu2sv3Adnp26nnYXwB1fRnomoEkOoW8JKTyinLyD+ZXhX7VF0CtxxVewaDu\ngzi2x7Fk9shkaK+hZPbMrNr6d+tPiumGNWm7FPKS1PYV7yN3Xy6bCzaTszfnsG3j3o0UFBWEL4BY\n6A/tefiXQHq3dH0JSFxTyIscwaHSQ2zau+nwL4CC6uN9xfuqvwR6ZHJcr+MY3nt41da1Q9eo/xMk\nySnkRZrgYMlBNhVUfwm8v+d9NuzewIbdG3hvz3v06tSL4b2HM6L3iKrgH9FnBMN6DSOtY1rUzZck\noJAXaSHP6hScAAAKRUlEQVQVXsHWfVvZsHsD63evrwr/9bvX897u9+jesTsj+oz40JfA8N7D6d4x\nWVa7kZamkBeJQIVXsG3/turg37WeDXs2VD3u1qHbh/8CiB336NQj6uZLG6KQF4kz7s4HBz6oDv/d\nG9iwp/o4rWMao44Zxci+Iw/b9+vaT8NB5UMU8iJtiLuTuy+X1TtWs2bnmsP27s7IY0Yyqu+osI99\nAQzuMVgjgJKYQl4kAbg7Ow7tYM2O6uCvDP+CogJO6HtCVehXfgEc1+s42qVo/d9Ep5AXSXAFRQUh\n9Gt9AWzbv43hvYcfVvYZecxIju9zPJ3adYq62dJMFPIiSepQ6SHW7VoXgn/HGlbvDPuNezcytOdQ\nxqSP4aR+J3FS+kmMSR/DkB5DVPNvgxTyInKYkvIS1u5cy8q8lazIW8GK/BWsyFvBgZIDjOk3hjH9\nxnBSegj/E/udqNE+ca7VQt7M7gI+AZQA7wHXu3tB7LkZwJeAcuAb7v5iHe9XyItEaNehXazMD8G/\nMm8lK/JX8E7+O/Tt0jf09muE/4g+I1TvjxOtGfJTgZfcvcLMbgdw9+lmNgp4BBgPZADzgOPdvaLW\n+xXyInGmvKKc9/e8XxX+K/JWsDJ/Jdv2b+PEficytv9YxvYfy7gB4xiTPka1/ghEUq4xs8uBT7n7\nNbFefIW73xF77nlgprsvqvUehbxIG7G/eD9v573Nsg+WsXz7cpZ/sJx1u9YxrPcwxg0YVxX+p/Q/\nReWeFtbQkG+uv7++BDwaOx4I1Az0XEKPXkTaqLSOaZx17FmcdexZVeeKy4pZlb+qKvRnvTOLlXkr\n6d+tP2MHjGVc/3GMHRDCP71beoStT25HDHkzmwv0r+OpH7r7M7HX/AgocfdHjvCj6uyyz5w5s+o4\nKyuLrKysozRXROJFx3YdOXXgqZw68NSqc+UV5azbta6qx3/XwrtY9sEyOrfrzLgB4zht4GmcNvA0\nxg8cr+Cvp+zsbLKzsxv9/iaVa8zsOuBG4Fx3L4qdmw7g7rfHHj8P3Obui2u9V+UakSTg7mwq2MSy\nD5bxxrY3qrZuHbpVBX5l+Pfq3Cvq5sa91rzweiFwN3C2u++scb7ywuvpVF94HV470RXyIsnL3Xl/\nz/ss3baUN7a9wdJtS1n2wTLSu6YzPmM8pw04jfEZ4xk3YBzdOnSLurlxpTVDfj3QAdgdO/W6u98S\ne+6HhDp9GfBNd3+hjvcr5EWkSnlFOe/uejeE/talLN22lJX5K8nsmcnpGaczIWMCEwdN5MR+Jyb1\ncE7dDCUiCaO0vJRV+atYsnUJi7cuZlHuIrbs28K4AeOqQn9CxgQyuifP2A6FvIgktL1Fe1m6dWlV\n6C/euphO7ToxIWNCVfCfOvBUurTvEnVTW4RCXkSSSmV9v2bor8pfxcf6fIyJgyZy5uAzmTR4Epk9\nMxNirh6FvIgkvaKyIpZ/sJxFuYtYmLuQBZsX4HhV4E8aPIlxA8bRIbVD1E1tMIW8iEgtlcM4F25Z\nyMItC1mwZQHrd61n7ICxTBo0qSr4j+l6TNRNPSqFvIhIPewv3s+SrUtYsGUBC7csZFHuIvp17cek\nwZOYfOxkpgyZwvDew+OuxKOQFxFphPKKctbsXMNrm1/j1c2v8sqmVygtL2XKkClMPnYyk4dMZky/\nMaSmpEbaToW8iEgzqCzxvLopBP6rm18l72AekwZPYsqxU5g8ZDKnDTyt1ev6CnkRkRaSdyDvsJ7+\nul3rGJ8xvqq8M3HQxBa/Q1chLyLSSgqKCli4ZWFV6C/fvpwT+51I1pAszhl6DmcOPpO0jmnN+pkK\neRGRiBSWFrJ462Kyc7KZnzOfN7e9yZj0MZyTeQ7nZJ7Dmcee2eSbtBTyIiJx4lDpIV7f8jrzc+Yz\nP2c+b29/m7EDxlaF/hmDz2jw6loKeRGROHWg5AALtyxk/sYQ+qvyVzE+YzznZJ5DVmYWEzIm0LFd\nxyP+DIW8iEgbsa94H69tfq0q9N/d9S4TMiZwTuY5fHzoxxmfMf5DM24q5EVE2qi9RXt5ddOrzM+Z\nz8sbXyZnbw6Th0zmvKHnce5x5zL6mNGkpKQo5EVEEkH+wXzmb5zPSxtfYt778ygsK2T797a3bsib\n2XeBu4C+7r47dm4GYdGQcuAb7v5iHe9TyIuINMD7e95nWO9hDQr5lKZ8oJkNBqYCm2qcGwVcBYwC\nLgTuNbMmfU6ia8oivYlGv4tq+l1U0+8iOK7XcQ1+T1PD91fAD2qdmwY86u6l7p4DbCCs9yofQf8D\nV9Pvopp+F9X0u2i8Roe8mU0Dct19Ra2nBgK5NR7nEhb0FhGRVnbE1XDNbC7Qv46nfgTMAM6v+fIj\n/CgV30VEItCoC69mdiLwEnAodmoQsBWYAFwP4O63x177PHCbuy+u9TMU/CIijdDqQyjNbCNwqrvv\njl14fYRQh88A5gHDNZRGRKT1HbFc0wBVAe7uq81sNrAaKANuUcCLiEQjspuhRESk5UUyft3MLjSz\ntWa23sxujaIN8cDMBpvZfDN7x8xWmdk3om5T1Mws1cyWm9kzUbclSmbW08weN7M1ZrbazCZG3aao\nmNmM2L+RlWb2iJkdeQavBGJm95tZnpmtrHGut5nNNbN1ZvaimfU80s9o9ZA3s1Tgd4QbpUYBnzOz\nka3djjhRCnzb3UcDE4GvJfHvotI3CaW+ZP8T89fAs+4+EjgJWBNxeyJhZpnAjcA4dx8DpAKfjbJN\nrewBQlbWNB2Y6+7HEwbATD/SD4iiJ386sMHdc9y9FHiMcANV0nH37e7+Vuz4AOEf8sBoWxUdMxsE\nXAz8mSMPyU1oZtYDmOzu9wO4e5m7F0TcrKjsI3SGuphZO6ALYSRfUnD3V4E9tU5fCjwYO34QuOxI\nPyOKkM8AttR4rJulqOqxjAUWH/mVCe0e4PtARdQNidhQYIeZPWBmy8zsPjNr2nJCbVRsPqy7gc3A\nNmCvu8+LtlWRS3f3vNhxHpB+pBdHEfLJ/mf4h5hZN+Bx4JuxHn3SMbNPAPnuvpwk7sXHtAPGAfe6\n+zjgIEf5kzxRmdkw4FtAJuGv3G5m9vlIGxVHYiMXj5ipUYT8VmBwjceDOXwahKRiZu2BJ4CH3f2p\nqNsToUnApbF7Lh4FPm5mD0XcpqjkEqYMWRp7/Dgh9JPRacBCd9/l7mXAHML/K8ksz8z6A5jZACD/\nSC+OIuTfAEaYWaaZdSDMWPmPCNoROTMz4C/Aanf/n6jbEyV3/6G7D3b3oYQLay+7+xejblcU3H07\nsMXMjo+dOg94J8ImRWktMNHMOsf+vZxHuDCfzP4BXBs7vhY4YuewuW6Gqjd3LzOzrwMvEK6U/8Xd\nk3LkAHAmcA2wwsyWx87NcPfnI2xTvEj2st5/AH+LdYTeIzZdSLJx97djf9G9QbhWswz4U7Staj1m\n9ihwNtDXzLYAPwFuB2ab2Q1ADnDlEX+GboYSEUlcWsxDRCSBKeRFRBKYQl5EJIEp5EVEEphCXkQk\ngSnkRUQSmEJeRCSBKeRFRBLY/wcdFqYC5H7KdgAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x6b59b10>"
       ]
      }
     ],
     "prompt_number": 97
    }
   ],
   "metadata": {}
  }
 ]
}