# Collision penetration problem

I'm implementing a 3D engine but I am having some problems with the physics. Sadly my education on this area is not that extend and i would love some help to solve this problem.

Let's say we have a infinite plane and a oriented box. The box has both linear and angular velocities, the plane is static, which means its linear and angular velocities are zero. In some moment, a collision happens, and we know which border, face or vertex of the cube collided with the plane and how depth it is. Let's illustrate it.

In this point we know a few things, let's list them:

• Angular velocity of the box in this frame (variation of rotation in this frame) (the axis of rotation is at the center of the box)
• Linear velocity of the center object of the box in this frame (variation of position in this frame)
• depth (it is a vector which direction is the normal of the plane and the magnitude is the vertex/border/face penetration)
• The position of the border/vertex/face that collided first with the plane after the collision (i called it collision point at the image)
• The normals of the box and the plane
• The position of center of mass of the box
• The inertia tensor of the box
• The mass of the box

I am looking for the time t in which the collision happened to get which is the angular velocity and linear velocity that the object had at that point.

Thank you very much for your attention :)

• Is this really a question about physics? Sep 29, 2016 at 10:22
• Do you know how linear and angular velocity change in time, or can they be assumed constant during this time interval? Sep 29, 2016 at 10:56
• They are constant. The angular and linear velocity are the variance of the position and rotation in one frame in this case. Sep 29, 2016 at 12:35
• It is a question about physics yes. How do you expect that a engine simulates reality if it doesn't use the same rules than reality? Sep 29, 2016 at 12:38
• @fibonatic the angular and linear velocities are the variance of those between t=0 and t=1. For that reason im looking for the t, because if i know the t i can multiple it with the linear and angular velocities to know the exact velocities of the objects when the collision happens. Sep 29, 2016 at 14:02

To find the intermediate position between two time steps you need to use a cubic spline function without any loss of precision. Let's say each body has three degrees of freedom $(x,y,\theta)$ and you know the values and their derivatives between two time steps $$\begin{pmatrix} x_1 \\ y_1 \\ \theta_1 \end{pmatrix},\begin{pmatrix} \dot{x}_1 \\ \dot{y}_1 \\ \dot{\theta}_1 \end{pmatrix} \rightarrow \begin{pmatrix} x_2 \\ y_2 \\ \theta_2 \end{pmatrix},\begin{pmatrix} \dot{x}_2 \\ \dot{y}_2 \\ \dot{\theta}_2 \end{pmatrix}$$

If the time step size is $h$ and $t=0$ is at the first step then the intermediate values are found by $$\begin{pmatrix} x \\ y \\ \theta \end{pmatrix} = \begin{vmatrix} x_1 & x_2 & \dot{x}_1 & \dot{x}_2 \\ y_1 & y_2 & \dot{y}_1 & \dot{y}_2 \\ \theta_1 & \theta_2 & \dot{\theta}_1 & \dot{\theta}_2 \\ \end{vmatrix} \begin{pmatrix} 2 \left( \tfrac{t}{h} \right)^3 - 3 \left( \tfrac{t}{h} \right)^2 + 1 \\ \left( \tfrac{t}{h} \right)^2 \left(3-2 \left( \tfrac{t}{h} \right) \right) \\ t \left( \left( \tfrac{t}{h} \right)^2 - 2 \left( \tfrac{t}{h} \right) +1 \right) \\ h \left( \tfrac{t}{h} \right)^2 \left(\tfrac{t}{h}-1\right) \end{pmatrix}$$

For each $(x,y,\theta)$ solution there are the coordinates of the contact corner P defined as $$(x_p,y_p) = f(x,y,\theta,\mbox{dimensions})$$

For a box of size $a$ and $b$ this is

$$\begin{pmatrix} x_p \\ y_p \end{pmatrix}= \begin{pmatrix} x + \frac{a}{2}\cos\theta-\frac{b}{2}\sin\theta \\ y + \frac{a}{2}\sin\theta + \frac{b}{2} \cos\theta \end{pmatrix}$$

The constraint is if the form of a line (plane in 3D) $n_x x + n_y y = d$ where $(n_x,n_y)$ is the contact normal vector and $d$ is the distance from the origin. You are at the contact point (within tolerance $\epsilon$) when

$$|n_x x_p + n_y y_p - d | \leq \epsilon$$

You can choose now to use a numerical method (like bisection) to solve the problem, or to assume the time step is small enough to solve it using the assumption that $t\ll 1$ and

$$\boxed{ t \approx \frac{d - n_x x_p -n_y y_p }{n_x (\dot{x}_1 +\dot{\theta}_1 (y_1 -y_p))-n_y ((x_1-x_p) \dot{\theta}_1-\dot{y}_1)} }$$

• It solved the problem! Thank you very much for your help :D. Oct 3, 2016 at 22:10

The fast and somewhat accurate way (depending on the size of the time step) to do this is to take the position just before the collision ($\vec{a}$) and just after the collision ($\vec{b}$), draw a line between them $$\vec{L}(t) = (1-t)\vec{a} + t\vec{b}$$ and solve for when $\vec{L}$ is in the plane. Namely, when $$(\vec{L}(t) - \vec{p})\cdot\vec{n} = 0,$$ where $\vec{p}$ is a point in the plane and $\vec{n}$ is the normal vector of the plane (if that's how you define planes).

For the physicsy solution, I'm only going to track a vertex since it's much easier to track a single point instead of a higher dimensional entity. Plus, an object that moving and rotating is much more likely to hit on a corner since an edge or face would require a precise alignment and timing.

More assumptions to make my life easier (mostly safe assumptions assuming $t = 0 \rightarrow 1$ is one frame of animation):

• Linear velocity is constant.
• Angular velocity is constant.
• Axis of rotation is constant.

We can describe the motion of the vertex as the combination of the motion of the center of mass and the rotational motion about the center of mass: $$\vec{x}(t) = \vec{x}_\textrm{CM}(t) + \vec{x}_\textrm R(t).$$ Center of mass is easy: $$\vec{x}_\textrm{CM}(t) = \vec{x}_\textrm{CM}(0) + \vec{v}t$$ where $\vec{v}$ is the velocity of the center of mass.

The change in position around the center of mass is a rotation of the vector $$\vec{x}' = \vec{x}(0) - \vec{x}_\textrm{CM}(0)$$ by an angle $\theta = \omega{}t$ about an axis vector $\vec{b}$, where $\omega$ is the angular velocity. The most straightforward way I've found to represent rotations about an arbitrary axis is a rotation matrix following a coordinate transform. The x-axis will be the unit vector in the direction of the component of $\vec{x}'$ perpendicular to the rotation axis $\vec{b}$ (call it $\hat{x}'_b$), the z-axis will be the unit vector of the rotation axis ($\hat{b}$), and the y-axis will be the cross product of those two. The inverse (and transpose) of this transform is given by $$D^{-1} = \begin{bmatrix} | & | & | \\ \hat{x}'_b & \hat{b} \times \hat{x}'_b & \hat{b} \\ | & | & | \end{bmatrix}$$ This makes the rotation a simple $$R = \begin{bmatrix} \cos(\theta) & -\sin(\theta) & 0 \\ \sin(\theta) & \phantom{-}\cos(\theta) & 0 \\ 0 & 0 & 1 \end{bmatrix}$$

The full transform is giving by \begin{align}D^{-1}RD\vec{x}' &= D^{-1}R\begin{bmatrix}|\vec{x}'_b|\\0\\\hat{b}\cdot\vec{x}'\end{bmatrix} \\ &= D^{-1}\begin{bmatrix}|\vec{x}'_b|\cos(\theta)\\|\vec{x}'_b|\sin(\theta)\\\hat{b}\cdot\vec{x}'\end{bmatrix} \\ &= \vec{x}'_b\cos(\theta) + \sin(\theta)(\hat{b} \times \vec{x}'_b) +(\hat{b}\cdot\vec{x}')\hat{b} \end{align}

So, the total movement is given by \begin{align} \vec{x}(t) &= \vec{x}_{CM}(t) + \vec{x}_R(t) \\ &= \vec{x}_{CM}(0) + \vec{v}t + \vec{x}'_b\cos(\omega{}t) + \sin(\omega{}t)(\hat{b} \times \vec{x}'_b) +(\hat{b}\cdot\vec{x}')\hat{b} \end{align} (note the changes of $\hat{x}'_b$ to $\vec{x}'_b$ due to the multiplication with $|\vec{x}'_b|$). To get $t$, solve for when $\vec{x}(t)$ is in the plane as in the first paragraph: $$\left(\vec{x}(t) - \vec{p}\right)\cdot\vec{n} = 0.$$

In this form, the equation is intractable with the $\sin$ and $\cos$. You can make the small angle approximations $$\sin x \approx x$$ and $$\cos x \approx 1 - \frac{1}{2}x^2$$ to turn it into a solvable quadratic equation.

If you don't have access to the position just before impact ($t=0$ position/orientation), you can use $$\left(\left(\vec{x}(1) - \vec{p}\right)\cdot\hat{n}\right)\hat{n} = \vec{d}$$ to find it (note the use of the unit normal vector $\hat{n}$).

• In the absence of forces each point on a body travels in a curvilinear path (rotation + translation). So only for very small time steps the linear approximation would work correctly. Sep 30, 2016 at 16:57
• @ja72 Right, which is why the answer continues after the first paragraph. The "physicsy" solution is the curvilinear path. Sep 30, 2016 at 19:58
• Thank you very much for your help! It solved the problem, but the solution made by @ja72 was faster for the computer to run, i would wish i could point both as correct :( Oct 3, 2016 at 22:08
• @Haruko I'm actually kinda surprised my answer worked, seeing as I was just typing thing from memory of other projects. Oct 6, 2016 at 10:50