Modelling a basketball 3 point swish I am trying to do some simulations of a 3 point basketball swish shot. A 'swish' is when the ball enters the hoop without touching the rim or backboard. I am confident with python but do not know how to proceed with this.
I would like to find the combinations of initial velocities $\vec v_0 = (v_{0x}, v_{0y}, v_{0z})$ and spin $\omega$ that allow the shooter to land a swish 3 pointer.  The position $(0,0,0)$ where the ball is shot. What I've done: -derived the equations of motion:
The velocity vector: $$ \vec v \begin{cases}
\dot{x}(t)=v_{x0}+\frac{F_M}{m}t-\frac{F_D}{m}\cdot\frac{\dot{x}(t)}{\sqrt{\dot{x}^2(t)+\dot{y}^2(t)+\dot{z}^2(t)}}t\\
\dot{y}(t)=v_{y0}-\frac{F_D}{m}\cdot\frac{\dot{y}(t)}{\sqrt{\dot{x}^2(t)+\dot{y}^2(t)+\dot{z}^2(t)}}t\\
\dot{z}(t)=v_{z0}-gt-\frac{F_D}{m}\cdot\frac{\dot{z}(t)}{\sqrt{\dot{x}^2(t)+\dot{y}^2(t)+\dot{z}^2(t)}}t
\end{cases}$$
The position vector: $$ \vec r = \begin{cases}
x=v_{x0}t+\frac{F_M t^2}{2m}-\frac{F_D t^2}{2m}\cdot\frac{\dot{x}(t)}{\sqrt{\dot{x}^2(t)+\dot{y}^2(t)+\dot{z}^2(t)}}\\
y=v_{y0}t-\frac{F_D t^2}{2m}\cdot\frac{\dot{y}(t)}{\sqrt{\dot{x}^2(t)+\dot{y}^2(t)+\dot{z}^2(t)}}\\
z=v_{z0}t-\frac{gt^2}{2}-\frac{F_D t^2}{2m}\cdot\frac{\dot{z}(t)}{\sqrt{\dot{x}^2(t)+\dot{y}^2(t)+\dot{z}^2(t)}} \end{cases} $$
-Derived some boundary conditions: I am assuming that at a time $T$ the centre of mass of the ball is at the same height as the rim so $z(T) = H$ and the other condition is that there needs to be a set of points $(x,y)$ which ensure that the sphere is enclosed within the radius of the rim and is not touching it. This constraint has the form $x^2(T) + (a-y(T))^2 \leq (r_{\mbox{ ring}} -r_b)^2$.
I would like to simulate this experiment and do some plots to determine the the combinations which allow me to land this shot. I am pretty solid with Python and have basic knowledge of Matlab. Any help would be amazing.
 A: Hint.
With $p=(x(t),y(t),z(t)),\ \ v = (v_x(t),v_y(t),v_z(t))$ assuming a dynamic model
$$
\cases{
\dot v = (0,0,-g)-\frac{F_D}{m}\frac{v}{||v||}\\
\dot p(t) = v(t)
}
$$
with initial conditions
$$
p(0) = (0,0,0),\ \ \  v(0) = (v_{x_0},v_{y_0},v_{z_0})
$$
Suppose we need to determine the orbit from $p_0=(0,0,0)$ to $p_{goal}=(a,b,H)$ with the restrictions $v_{x_0}^2+v_{y_0}^2+v_{z_0}^2\le E_0$ and $\frac{\dot p(T)}{||\dot p(T)||}\cdot (0,0,-1)\le \cos\left(\alpha_0\right)$. Our approach is to handle this as a minimization problem with the objective function
$$
f=\left(\frac{v(T)}{||v(T)||}\cdot (0,0,-1)-\cos\left(\alpha_0\right)\right)^2+\left(p(T)-p_{goal}\right)\cdot \left(p(T)-p_{goal}\right)
$$
Follows a MATHEMATICA script which implements this.
Clear[Goal, vx, vy, vz, T, p, v]
p = {x[t], y[t], z[t]};
v = {vx[t], vy[t], vz[t]};
pgoal = {a, a, H};
modv = Sqrt[v.v];
parms = {Fd -> 1/10, m -> 1, g -> 10, a -> 10, E0 -> 180, H -> 2.2, alpha0 -> Pi/5};

Goal[vx0_?NumberQ, vy0_?NumberQ, vz0_?NumberQ, T_?NumberQ] := ((v/modv).{0, 0,-1} - Cos[Pi/5])^2 + (p - pgoal).(p - pgoal) /. 
NDSolve[{vx'[t] == -Fd/m vx[t]/modv, 
         vy'[t] == -Fd/m vy[t]/modv, 
         vz'[t] == -g - Fd/m vz[t]/modv, 
         x'[t] == vx[t], 
         y'[t] == vy[t], 
         z'[t] == vz[t], 
         vx[0] == vx0, vy[0] == vy0, vz[0] == vz0, 
         x[0] == 0, y[0] == 0, z[0] == 0} /. parms, {x, y, z, vx, vy, vz}, {t, 0, T}][[1]] /. parms /. {t -> T}

And a minimization procedure
sol = NMinimize[{Goal[vx0, vy0, vz0, T], vx0^2 + vy0^2 + vz0^2 <= E0,T > 0}, {vx0, vy0, vz0, T}, Method -> "DifferentialEvolution"]

And the graphical verification
sol0 = NDSolve[{vx'[t] == -Fd/m vx[t]/modv, 
                vy'[t] == -Fd/m vy[t]/modv, 
                vz'[t] == -g - Fd/m vz[t]/modv, 
                x'[t] == vx[t], 
                y'[t] == vy[t], 
                z'[t] == vz[t], 
                vx[0] == vx0, vy[0] == vy0, vz[0] == vz0, 
                x[0] == 0, y[0] == 0, z[0] == 0} /. parms /. sol[[2]], {x, y, z, vx, vy, vz}, {t, 0, T /. sol[[2]]}][[1]]

ParametricPlot3D[Evaluate[p /. sol0], {t, 0, T /. sol[[2]]}]


