__author__ = 'vladimirtsvetkov'
import numpy as np

a = np.pi / 4
g = 10
v0 = 40
ro = 7600
r = 0.2

#first model(simple)
x = np.tan(a) / (g / (2 * v0 * v0 * np.cos(a) * np.cos(a)))
print 'simple model:'
print x

#second model(complex)
vx = v0 * np.cos(a)
vy = v0 * np.sin(a)
h = 0.01
b = 0.16 * 1.15 * np.pi * r * r / 2.0
m = ro * 4 * np.pi * r * r * r /3

def fvx(vx, vy):
    return ((-b * vx * np.sqrt(vx*vx + vy*vy))/m)

def fvy(vx, vy):
    return -g - ((b * vy * np.sqrt(vx*vx + vy*vy))/m)

y = 0
x = 0
flag = False
xprev = 0
yprev = 0

xmax = 120

while (True):
    k11 = h * fvx(vx, vy)
    k12 = h * fvy(vx, vy)

    k21 = h * fvx(vx + h * k11/2, vy + h * k12/2)
    k22 = h * fvy(vx + h * k11/2, vy + h * k12/2)

    k31 = h * fvx(vx + h * k21/2, vy + h * k22/2)
    k32 = h * fvy(vx + h * k21/2, vy + h * k22/2)

    k41 = h * fvx(vx + h * k31, vy + h * k32)
    k42 = h * fvy(vx + h * k31, vy + h * k32)

    vx += (k11 + 2 * k21 + 2 * k31 + k41)/6
    vy += (k12 + 2 * k22 + 2 * k32 + k42)/6

    v = np.sqrt(vx*vx + vy*vy)

    xprev = x
    yprev = y

    x += vx * h
    y += vy * h

    if (x == xmax):
        flag = True
        break


    if (x > xmax):
        break

if (flag):
    print 'complex model1:'
    print x
else:
    print 'complex model2:'
    print x
    #print (xprev * y - x * yprev)/(y - yprev)