How to correct for collision in $N$-dimensional bounded motion of a particle? [closed]

Question

I want to numerically simulate the motion of a particle (unit mass) in $N$-dimensional space according to the force $f(\vec{x})$ and bounded by a set of linear constraints (equalities as well as inequalities):

$$A_i \vec{x} = \vec{c_i}$$

$$B_j \vec{x} \leq \vec{d_j}$$

The behavior of the particle when it comes in contact with a "wall" formed by an inequality is a perfectly inelastic collision.

What I'm thinking right now is something like this:

for step in range(total_time//dt):
    p += f(x) * dt
    x += p * dt
    if violate_constraints(x):
        x, p = correct(x, p)

It is easy to check if $\vec{x}$ violates any constraints, but I'm not sure how to implement the correct(x, p) part.

I think the correct behavior for correct(x, p) is to put the position $\vec{x}$ back to the closest point in the "allowed" space, and project $\vec{p}$ onto the intersection of subspaces according to which constraints $\vec{x}$, before correction, broke. However,

1. I have no proof if this is correct
2. I don't know how to put $\vec{x}$ back to the closet point in the allowed space

EDIT: Putting $\vec{x}$ back to the closet point in the allowed space is what I think should correctly handle the "corner" cases where the particle comes in contact with 2 or more walls.

Answer 1

Let's say you write your vectors as lists with length n. Then I would write something like this

dt=...
x=[...] #your initial state
p=[...]
def f(x):
    f = []
    f.append(${f_i(x)})
    return f

def evolve_to_next_state(x, p):
    previous_x = x
    previous_p = p
    for d in range(len(x)):
        p[d] += f(x)[d]*dt
        x[d] += p[d]*dt
    for d in range(len(x)):
        if violate_constraints(x[d]):
            x[d], p[d] = correct(d, previous_x, previous_p)

 def correct(d, previous_x, previous_p):
    x[d] = previous_x[d]
    p[d] = -previous_p[d]
    return x[d], p[d]

This should solve the corner problem, going through each dimension singularly. Don't forget that now violate_constraints takes floats as input. You then have to write a function that plots it, or at least goes through each time interval... Hope it could help

Edit: 1. assume point particle... 2. Obviously this corrects only the inequality, you have to implement the code to correct also the linear equality condition

Answer 2

With a perfectly inelastic collision between the wall and the mass, the mass will stick to the wall at the point of the collision and all of its velocity drops to 0. That part of the problem should be straight forward to implement (e.g., if collision_detected: v = 0). The goal is then to find the point at which it collided with the wall and set its position to that point.

Using the image below, produced by Nikos M on this StackOverflow post, we can see that if $x_\text{mass}>x_\text{wall}$, then we can compute the depth, $D$, past the wall's edge via, $$ D=x_\text{mass}-x_\text{wall}+r. $$ where I assume the mass is a ball of radius $r$ (obviously if it's a point mass, then $r\to0$ is simple enough a task).

You still have your previous position increment, $p\cdot\mathrm{d}t$, so your $x$ distance is reduced from $p_x\cdot\mathrm{d}t$ to $p_x\cdot\mathrm{d}t-D$. Simple geometry should lead you to the angle of attack, $\tan\theta=p_y/p_x$, so then you can determine the height you can go, $$H=p_y\cdot\mathrm{d}t-p_x\cdot\mathrm{d}t\tan\theta$$ Thus, you have your new positions when crossing the wall along the $x$ axis: $\vec{x}_\text{new}=D\hat{\mathrm{e}}_x+H\hat{\mathrm{e}}_y$.

(source)

The case for crossing the $y$ axis is done with the same procedure. For the corner cases, it seems to me that you can simply place the mass at the corner itself, that is, $$\vec{x}_\text{new}=(x_\text{wall}-r)\hat{\mathrm{e}}_x+(y_\text{wall}-r)\hat{\mathrm{e}}_y,$$ as the collision with the two walls would place it directly in the corner, rather than trying to do a sweep in the $x$ and then $y$ direction (though that might be necessary for cases where the collision is near-enough to the corner, I'm not sure).