I'm trying to create an eight-ball pool simulation in Python, and I'm trying to simulate the interaction between a cue stick and a cue ball. Suppose I thrust a cue stick with mass $m_s$ at a velocity $v_s$ to hit a cue ball at a position $(a, b)$ from its center in the $xz$ plane. I also know the ball's mass, $m_b$, and its radius, $R$. I would like to determine the velocity components $\{v_x, v_y, v_z\}$ and angular velocity components $\{\omega_x, \omega_y, \omega_z\}$ after the cue stick hits the cue ball.
Here is a picture describing the interaction between the cue stick and the ball in the above example:
I have come up with some equations using the help of this article, but my results seem incorrect after comparing them with the results of another pool simulator. I would like to confirm whether the equations (highlighted in bold at the bottom) are valid, and if not, adjust them such that they are valid.
By the conservation of linear momentum: $$m_sv_{s_i} = m_sv_{s_f} + m_bv_b \tag{1}\label{1}$$
Assuming an elastic collision, by the conservation of energy (the kinetic energy of the cue is converted into the final kinetic energy of the cue, the kinetic energy of the cue ball, and the rotational kinetic energy of the cue ball): $$\frac{1}{2}m_s(v_{s_i})^2 = \frac{1}{2}m_s(v_{s_f})^2 + \frac{1}{2}m_b(v_b)^2 + \frac{1}{2}(\frac{2}{5}m_bR^2)(\omega_x)^2 + \frac{1}{2}(\frac{2}{5}m_bR^2)(\omega_z)^2 \tag{2}\label{2}$$
Since the linear impulse between the tip and the ball is equal to the change in momentum of both the stick and the ball: $$Ft = m_s(v_{s_i} - v_{s_f}) = m_bv_b \tag{3}\label{3}$$
Since the angular impulse is equal to the change in angular momentum: $$bFt = (\frac{2}{5}m_bR^2)(\omega_x) \tag{4}\label{4}$$ $$aFt = (\frac{2}{5}m_bR^2)(\omega_z) \tag{5}\label{5}$$
Then, using $(1)$, the final stick speed is: $$v_{s_i} - \frac{m_b}{m_s}v_b = v_{s_f} \tag{6}\label{6}$$
Then, substituting $(3)$ into $(4)$ and $(3)$ into $(5)$, $$bm_bv_b = (\frac{2}{5}m_bR^2)(\omega_x) \\ \frac{5}{2}\frac{v_b}{R^2}b = \omega_x \tag{7}\label{7}$$ and similarly:
$$\frac{5}{2}\frac{v_b}{R^2}a = \omega_z \tag{8}\label{8}$$
Substituting $(6)$, $(7)$, and $(8)$ into (2) yields: $$v_b = \frac{2v_s}{1 + \frac{m_b}{m_s} + \frac{5}{2}(\frac{a+b}{R})^2} \tag{9}\label{9}$$
Since $v_b\sin(\theta) = v_x$ and $v_b\cos(\theta) = v_y$, we have the velocity components $\{v_b\sin(\theta), v_b\cos(\theta), 0\}$ and the angular velocity components $\{\frac{5}{2}\frac{v_b}{R^2}b, 0, \frac{5}{2}\frac{v_b}{R^2}a\}$, where $v_b$ is obtained from equation $(9)$.
EDIT:
Dug a little deeper and found what I was looking for in a paper titled Pool Physics Simulation by Event Prediction I. In this article, the authors consider the cue angled to the horizontal $xy$ plane. Assuming the time duration of the collision between the cue stick and the ball is small, the force exerted on the ball can be "treated as a perfectly elastic impulse" and thus the velocity of the ball can be expressed as: $$v= (0, -\frac{F}{m}\cos(\theta), -\frac{F}{m}\sin(\theta))$$
The authors claim that "the magnitude of the force $F$ in terms of the impact parameters is derived simultaneously solving the equations for the conservation of linear momentum and conservation of energy before and after the cue impact to obtain:" $$F=\frac{2mV_0}{1+\frac{m}{M}+\frac{5}{2R^2}(a^2 + b^2\cos^2(\theta) + c^2\sin^2(\theta) -2bc\cos(\theta)\sin(\theta))}$$
where $m$ is the ball mass and $M$ is the cue mass and $c = |\sqrt{R^2 - a^2 - b^2}|$.
How did the authors derive this equation for force?
Furthermore, this pool simulation claims that the numerator should be $2MV_0$. Which one is correct?
Thanks in advance.