{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import k3d\n",
    "import numpy as np\n",
    "import time\n",
    "\n",
    "from scipy import integrate\n",
    "\n",
    "plot = k3d.plot()\n",
    "\n",
    "line = k3d.line([[0,0,0]])\n",
    "\n",
    "plot += line"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ipywidgets import interact, interactive, fixed\n",
    "import ipywidgets as widgets\n",
    "\n",
    "\n",
    "@interact(Cx=widgets.FloatSlider(value=0, min=0, max=2.0))\n",
    "def g(Cx):\n",
    "    wiatr_x = -21.1\n",
    "    wiatr_y = -12.1\n",
    "    wiatr_z = 0.01\n",
    "    n = 200\n",
    "    g = 9.81\n",
    "    x0, y0, z0, vx0, vy0,vz0 = [0,0,0,10,0,10]\n",
    "    dt = 2.1/n\n",
    "    trajektoria = [ (x0,y0,z0) ]\n",
    "    for i in range(n):\n",
    "        vx = vx0 - Cx*(vx0-wiatr_x)*dt\n",
    "        vy = vy0 - Cx*(vy0-wiatr_y)*dt\n",
    "        vz = vz0  - g * dt - Cx*(vz0-wiatr_z)*dt\n",
    "        x = x0 + vx0 *dt\n",
    "        y = y0 + vy0 *dt\n",
    "        z = z0 + vz0 *dt\n",
    "        if z<0:\n",
    "            break\n",
    "        x0, y0, z0, vx0, vy0, vz0 = x, y, z, vx, vy, vz\n",
    "        trajektoria.append(( x,y,z ))\n",
    "    line.vertices = trajektoria\n",
    "\n",
    "plot.display()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "plot.camera_auto_fit=False\n",
    "plot.grid_auto_fit = False"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Double pendulum\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "# %load http://matplotlib.org/examples/animation/double_pendulum_animated.py\n",
    "# Double pendulum formula translated from the C code at\n",
    "# http://www.physics.usyd.edu.au/~wheat/dpend_html/solve_dpend.c\n",
    "\n",
    "from numpy import sin, cos\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import scipy.integrate as integrate\n",
    "import matplotlib.animation as animation\n",
    "\n",
    "G = 9.8  # acceleration due to gravity, in m/s^2\n",
    "L1 = 1.0  # length of pendulum 1 in m\n",
    "L2 = 1.0  # length of pendulum 2 in m\n",
    "M1 = 1.0  # mass of pendulum 1 in kg\n",
    "M2 = 1.0  # mass of pendulum 2 in kg\n",
    "\n",
    "\n",
    "def derivs(state, t):\n",
    "\n",
    "    dydx = np.zeros_like(state)\n",
    "    dydx[0] = state[1]\n",
    "\n",
    "    del_ = state[2] - state[0]\n",
    "    den1 = (M1 + M2)*L1 - M2*L1*cos(del_)*cos(del_)\n",
    "    dydx[1] = (M2*L1*state[1]*state[1]*sin(del_)*cos(del_) +\n",
    "               M2*G*sin(state[2])*cos(del_) +\n",
    "               M2*L2*state[3]*state[3]*sin(del_) -\n",
    "               (M1 + M2)*G*sin(state[0]))/den1\n",
    "\n",
    "    dydx[2] = state[3]\n",
    "\n",
    "    den2 = (L2/L1)*den1\n",
    "    dydx[3] = (-M2*L2*state[3]*state[3]*sin(del_)*cos(del_) +\n",
    "               (M1 + M2)*G*sin(state[0])*cos(del_) -\n",
    "               (M1 + M2)*L1*state[1]*state[1]*sin(del_) -\n",
    "               (M1 + M2)*G*sin(state[2]))/den2\n",
    "\n",
    "    return dydx\n",
    "\n",
    "# create a time array from 0..100 sampled at 0.05 second steps\n",
    "dt = 0.015\n",
    "t = np.arange(0.0, 20, dt)\n",
    "\n",
    "# th1 and th2 are the initial angles (degrees)\n",
    "# w10 and w20 are the initial angular velocities (degrees per second)\n",
    "th1 = 120.0\n",
    "w1 = 0.0\n",
    "th2 = -10.0\n",
    "w2 = 0.0\n",
    "\n",
    "# initial state\n",
    "state = np.radians([th1, w1, th2, w2])\n",
    "\n",
    "# integrate your ODE using scipy.integrate.\n",
    "y = integrate.odeint(derivs, state, t)\n",
    "\n",
    "x1 = L1*sin(y[:, 0])\n",
    "y1 = -L1*cos(y[:, 0])\n",
    "\n",
    "x2 = L2*sin(y[:, 2]) + x1\n",
    "y2 = -L2*cos(y[:, 2]) + y1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import k3d\n",
    "plot = k3d.plot(antialias=True)\n",
    "plot.display()\n",
    "\n",
    "plot.grid_auto_fit = False\n",
    "plot.canera_auto_fit = False\n",
    "configuration = np.array([[[0,0,0]]])\n",
    "double_pendulum = k3d.line(configuration,color=0xFF0000 ,width=5)\n",
    "double_pendulum_masses = k3d.line(configuration,color=0x0000ff ,width=5)#point_size=10)\n",
    "\n",
    "plot +=  double_pendulum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "plot.camera_auto_fit=False\n",
    "import time\n",
    "for i in range(len(x1)):\n",
    "    X = np.array( [[0,0,2],[x1[i],0,y1[i]+2],[x2[i],0,y2[i]+2]] )\n",
    "    double_pendulum.vertices = X    \n",
    "    time.sleep(0.01)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from ipywidgets import interact, interactive, fixed\n",
    "import ipywidgets as widgets\n",
    "from IPython.display import clear_output\n",
    "\n",
    "@interact(i=widgets.IntSlider(value=0,min=0,max=x1.shape[0]-1))\n",
    "def g(i):\n",
    "    X = np.array( [[0,0,2],[x1[i],0,y1[i]+2],[x2[i],0,y2[i]+2]] )\n",
    "    double_pendulum.vertices = X[:]\n",
    "    clear_output(wait=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}