Projectile Motion with Differential Equations

Question

I've been trying to write code that will calculate the required intercept angles of a projectile launch to hit a moving target, a person running around a game world. I have no issues programming things, but I need a bit of help with the mathematics.

The projectile's motion in a 3-dimensional Cartesian coordinate system is modeled as the solutions to the following initial value problems (taking advantage of the fact that motion along one axis is independent of the others). The coordinate system is rotated so that y points "upward" (against the direction of gravity), and x and z point "horizontally". It is still right-handed. Sorry about the ASCII equations; I don't have enough reputation to post images.

$x(0) = c_1$
$x'(0) = v \cos(θ)\cos(φ)$
$x''(t) = -0.01x'(t)$

$y(0) = c_3$
$y'(0) = v \sin(φ)$
$y''(t) = -0.03 + -0.01y'(t)$

$z(0) = c_2$
$z'(0) = v \sin(θ)\cos(φ)$
$z''(t) = -0.01x'(t)$

θ (horizontal angle), φ (elevation angle), and v are the launch parameters, but I can only control the two angles.

I have no problem solving these differential equations to derive the components of the projectile's velocity, but I have a moving target to contend with. Let's say X(t, θ, φ, v) is the projectile's position as a function of time with given launch parameters guaranteed to stay constant with respect to changing t. I also have Y(t) as the target's position as a function of time, which is differentiable. However, this function doesn't have a clean analytical solution because the moving person's acceleration is often not differentiable (as when trying to intercept someone jumping off a cliff -- his acceleration will be about -0.03 m/s^2 until he hits the ground, when it suddenly goes to 0). I know I can remove this discontinuity, but still I'm left with a piecewise function which is really hard to use. This made my second approach by modeling the relative positions, velocities, and accelerations with differential equations fail.

I've also tried simplifying the x and z coordinates into just one coordinate, and modeling the projectile as traveling an arc within a plane; this works well since there's no wind, and I can hit any static target in range accurately. But this only calculates based on the target's position and essentially calculates a function Z(φ) that is the altitude of the projectile when it is either directly above or below the target. I then calculate the roots of Z by using Newton's method, getting the angle (sometimes two).

I tried adapting the static method by calculating the time-of-flight of the projectile to the target's original position, using this to estimate the target's moved position, repeating for a few iterations, and then firing at the final moved position, but this method resulted in worse accuracy than simply spreading shots out around the static target. I think it diverged.

Now I'm stuck. I can't figure out a numerical method to find the required launch parameters given the target's position Y(t) as a function of time and a function X(t, θ, φ, v) and their derivatives. Is there anything that gives a reasonable degree of accuracy (10% or so) without anything a modern Intel processor can't do in about 50ms?

Edit I don't care when the time of impact happens; I just need the paths to intersect at a given time. I can calculate the target's position at any point in the near future (or in the past).

Answer 1

I think the idea of combining $x$ and $z$ into a single lateral coordinate is a good one. So you have a plane whose coordinates I'll call $r$ and $h,$ where $r = \sqrt{(x - x(0))^2 + (z - z(0))^2}$ is the lateral distance from the point where you launch the projectile and $h = y - y(0)$ is the relative height above that point.

For any given initial speed $v,$ at every point in that plane there is some minimum time at which that point can be hit by a projectile. That is, this time is a function of $v,$ $r,$ and $h$: call it $\tau(v,r,h).$ The launch angle $\varphi$ that permits hitting that point at that time is also a function of $v,$ $r,$ and $h,$ $\varphi(v,r,h).$

You could put the target's motion onto the same two coordinates, that is, express its motion in terms of the two functions $r = r_T(t)$ and $h = h_T(t).$ Then $f(t) = \tau(v, r_T(t), h_T(t)) - t)$ is a function of $t$ alone. Find a zero of $f.$ Having found $t_1$ such that $f(t_1) = 0,$ the angle of elevation to hit the target is $\varphi(v,r_T(t_1),h_T(t_1))$ and the azimuth $\theta$ can be derived from the components of horizontal displacement at that time, $x_T(t_1) - x(0)$ and $z_T(t_1) - z(0).$

The trick to doing this quickly is to calculate $\tau(v,r,h)$ quickly. I suggest a large precomputed three-dimensional table giving values of $\tau(v,r,h)$ at many possible values of $v,$ $r,$ and $h.$ For the actual values of $v,$ $r,$ and $h$ that you encounter while trying to find a zero of $f,$ you interpolate within the table. If you make the "steps" between input coordinates in your table small enough, linear interpolation should be OK. The question of linear interpolation in an $n$-dimensional table ($n=3$ in this case) has been addressed in Formula for $N$-Dimensional linear interpolation and Estimating line paths in vector fields. (You can also search the Web for "multilinear interpolation".)

Answer 2

I'm a bit confused regarding your edit: "I don't care when the time of impact happens; I just need the paths to intersect at a given time." Once you have $\mathbf{Y}(t)$ and $\mathbf{X}(t,\theta,\phi,v)$, determining inputs $\theta$ and $\phi$ such that $\mathbf{Y}(t_{c}) = \mathbf{X}(t_{c},\theta,\phi,v)$ for some collision time $t_{c}$ is just algebraic manipulation and equation solving.

For example, your model has the following solution for $x(t)$: $$x(t) = c_{1} + 100v\cos(\theta)\cos(\phi)(1-e^{-.01t}).$$ Call the target's $x$ coordinate $x_{T}$. Given a desired collision time $t_{c}$, simply set $x(t_{c})=x_{T}(t_{c})$: $$x(t_{c}) = x_{T}(t_{c})\rightarrow \cos(\theta)\cos(\phi) = \frac{x_{T}(t_{c})-c_{1}}{100v(1-e^{-.01t_{c}})}.$$ Once you solve $y_{T}(t_{c}) = y(t_{c})$ and $z_{T}(t_{c}) = z(t_{c})$ similarly, then you'll have enough information to solve explicitly for $\phi$ and $\theta$, and you're done.

Note that since this is a control problem, you may also be experiencing a variety of issues. For starters, your ability to plan a collision is highly dependent on the trajectories and velocities of the projectile and target. Another thing to consider is that for any $t_{c}$, there is definitely a maximum threshold time $t^{*}<t_{c}$ past which planning a collision for $t_{c}$ is no longer feasible.

Iteration time of the collision planning algorithm is another crucial parameter that requires sophisticated tuning, because if the control algorithm operates too frequently then stability issues can arise, but if it is operated too infrequently you're more likely to miss your threshold time and be unable to make the collision happen.

Although the latter issue can only occur in more realistic modeling in which the simulation is still active while the control algorithm is operating. So you can probably ignore it. Please let me know if you have any other questions, or if I missed what you're asking for. I have a lot of experience with control-based collision modeling and I'd love to help.