def newton alpha beta if alpha or alpha pi print Bad angle return def

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
def newton(v, alpha, beta, m):
if alpha <= 0 or alpha >= pi/2:
print("Bad angle")
return .0
def fx(t, v):
return -(beta*v**2)/m
def fy(t, v):
return -copysign(1, v)*beta*v**2/m-g
du, dw = RK4(fx), RK4(fy)
dt = 0.0001
u, w = v*cos(alpha), v*sin(alpha)
t = dt
x_old = y_old = 0
x, y = u*dt, w*dt
while (y > 0):
t += dt
u_new, w_new = u + du(t, u, dt), w + dw(t, w, dt)
x_old = x
x += (u_new + u) * dt / 2
y_old = y
y += (w_new + w) * dt / 2
u, w = u_new, w_new
if y_old == 0:
return x_old
return (x_old*(y_old - y) - y_old*(x_old-x))/(y_old - y)