In [ ]:

    

import k3d
import numpy as np
import time

from scipy import integrate

plot = k3d.plot()

line = k3d.line([[0,0,0]])

plot += line


    

In [ ]:

    

from ipywidgets import interact, interactive, fixed
import ipywidgets as widgets


@interact(Cx=widgets.FloatSlider(value=0, min=0, max=2.0))
def g(Cx):
    wiatr_x = -21.1
    wiatr_y = -12.1
    wiatr_z = 0.01
    n = 200
    g = 9.81
    x0, y0, z0, vx0, vy0,vz0 = [0,0,0,10,0,10]
    dt = 2.1/n
    trajektoria = [ (x0,y0,z0) ]
    for i in range(n):
        vx = vx0 - Cx*(vx0-wiatr_x)*dt
        vy = vy0 - Cx*(vy0-wiatr_y)*dt
        vz = vz0  - g * dt - Cx*(vz0-wiatr_z)*dt
        x = x0 + vx0 *dt
        y = y0 + vy0 *dt
        z = z0 + vz0 *dt
        if z<0:
            break
        x0, y0, z0, vx0, vy0, vz0 = x, y, z, vx, vy, vz
        trajektoria.append(( x,y,z ))
    line.vertices = trajektoria

plot.display()


    

In [ ]:

    

plot.camera_auto_fit=False
plot.grid_auto_fit = False


    

In [ ]:

    

%matplotlib inline

# %load http://matplotlib.org/examples/animation/double_pendulum_animated.py
# Double pendulum formula translated from the C code at
# http://www.physics.usyd.edu.au/~wheat/dpend_html/solve_dpend.c

from numpy import sin, cos
import numpy as np
import matplotlib.pyplot as plt
import scipy.integrate as integrate
import matplotlib.animation as animation

G = 9.8  # acceleration due to gravity, in m/s^2
L1 = 1.0  # length of pendulum 1 in m
L2 = 1.0  # length of pendulum 2 in m
M1 = 1.0  # mass of pendulum 1 in kg
M2 = 1.0  # mass of pendulum 2 in kg


def derivs(state, t):

    dydx = np.zeros_like(state)
    dydx[0] = state[1]

    del_ = state[2] - state[0]
    den1 = (M1 + M2)*L1 - M2*L1*cos(del_)*cos(del_)
    dydx[1] = (M2*L1*state[1]*state[1]*sin(del_)*cos(del_) +
               M2*G*sin(state[2])*cos(del_) +
               M2*L2*state[3]*state[3]*sin(del_) -
               (M1 + M2)*G*sin(state[0]))/den1

    dydx[2] = state[3]

    den2 = (L2/L1)*den1
    dydx[3] = (-M2*L2*state[3]*state[3]*sin(del_)*cos(del_) +
               (M1 + M2)*G*sin(state[0])*cos(del_) -
               (M1 + M2)*L1*state[1]*state[1]*sin(del_) -
               (M1 + M2)*G*sin(state[2]))/den2

    return dydx

# create a time array from 0..100 sampled at 0.05 second steps
dt = 0.015
t = np.arange(0.0, 20, dt)

# th1 and th2 are the initial angles (degrees)
# w10 and w20 are the initial angular velocities (degrees per second)
th1 = 120.0
w1 = 0.0
th2 = -10.0
w2 = 0.0

# initial state
state = np.radians([th1, w1, th2, w2])

# integrate your ODE using scipy.integrate.
y = integrate.odeint(derivs, state, t)

x1 = L1*sin(y[:, 0])
y1 = -L1*cos(y[:, 0])

x2 = L2*sin(y[:, 2]) + x1
y2 = -L2*cos(y[:, 2]) + y1


    

In [ ]:

    

import k3d
plot = k3d.plot(antialias=True)
plot.display()

plot.grid_auto_fit = False
plot.canera_auto_fit = False
configuration = np.array([[[0,0,0]]])
double_pendulum = k3d.line(configuration,color=0xFF0000 ,width=5)
double_pendulum_masses = k3d.line(configuration,color=0x0000ff ,width=5)#point_size=10)

plot +=  double_pendulum


    

In [ ]:

    

plot.camera_auto_fit=False
import time
for i in range(len(x1)):
    X = np.array( [[0,0,2],[x1[i],0,y1[i]+2],[x2[i],0,y2[i]+2]] )
    double_pendulum.vertices = X    
    time.sleep(0.01)


    

In [ ]:

    

from ipywidgets import interact, interactive, fixed
import ipywidgets as widgets
from IPython.display import clear_output

@interact(i=widgets.IntSlider(value=0,min=0,max=x1.shape[0]-1))
def g(i):
    X = np.array( [[0,0,2],[x1[i],0,y1[i]+2],[x2[i],0,y2[i]+2]] )
    double_pendulum.vertices = X[:]
    clear_output(wait=True)


    

In [ ]: