Compare commits

..

10 Commits

Author SHA1 Message Date
Karma Riuk
b36382d6ec Added final pdfs 2023-08-02 12:08:49 +02:00
Karma Riuk
567ac11dc6 applyed edits from Prof. Carzaniga 2023-08-02 12:08:41 +02:00
Karma Riuk
b44c6e5fb0 fixed minor mistake in formula 2023-08-02 12:08:12 +02:00
Karma Riuk
af13b1fd5d Final compilation to pdf 2023-06-22 18:20:02 +02:00
Karma Riuk
75c272bc89 Completed the computation of the impulse parameter 2023-06-22 18:19:44 +02:00
Karma Riuk
95fa9cf725 Finished the conclusion 2023-06-22 18:19:31 +02:00
Karma Riuk
a0a5f2f624 Finished the implementation section 2023-06-22 18:19:12 +02:00
Karma Riuk
9879a04efe Added some labels 2023-06-22 18:18:49 +02:00
Karma Riuk
0b1274b92a Written the collision resolution section 2023-06-14 11:59:21 +02:00
Karma Riuk
a567c6571c Fixed figures that were old 2023-06-13 22:57:14 +02:00
15 changed files with 569 additions and 61 deletions

1
BP_Arnaud_Fauconnet.pdf Symbolic link
View File

@ -0,0 +1 @@
bachelorproject.pdf

Binary file not shown.

View File

@ -5,25 +5,45 @@
\usepackage[]{subcaption} \usepackage[]{subcaption}
\usepackage[]{float} \usepackage[]{float}
\usepackage[]{multicol} \usepackage[]{multicol}
\usepackage{mathtools}
\usepackage[]{listings}
\usepackage{tikz} \usepackage{tikz}
\usepackage{tkz-euclide} \usepackage{tkz-euclide}
\usetikzlibrary{external,shapes,through} \usetikzlibrary{external,shapes,through,arrows}
\tikzexternalize[prefix=figures/] \tikzexternalize[prefix=figures/]
\tikzstyle{none}=[] \tikzstyle{none}=[]
\input{./tikzs/styles.tikzstyles} \input{./tikzs/styles.tikzstyles}
\tikzset{>=stealth}
\pgfdeclarelayer{nodelayer} \pgfdeclarelayer{nodelayer}
\pgfdeclarelayer{edgelayer} \pgfdeclarelayer{edgelayer}
\pgfsetlayers{edgelayer,nodelayer,main} \pgfsetlayers{edgelayer,nodelayer,main}
\lstset{language=C++,
basicstyle=\ttfamily,
keywordstyle=\color{blue}\ttfamily,
stringstyle=\color{red}\ttfamily,
commentstyle=\color{green}\ttfamily,
morecomment=[l][\color{magenta}]{\#},
emph={std, vector, vec2d, polygon, collision},
emphstyle=\color{orange}\ttfamily,
numbers=left,
morekeywords={uint},
numberstyle=\footnotesize,
numbersep=5pt,
frame=lines,
% breaklines=true,
% breaklines=true,
% postbreak=\raisebox{0ex}[0ex][0ex]{\space\ensuremath{\hookrightarrow}},
}
\graphicspath{{../figures/}{./figures/}} \graphicspath{{../figures/}{./figures/}}
\newcommand*{\vv}[1]{\overrightarrow{#1}} \newcommand*{\vv}[1]{\overrightarrow{#1}}
\newcommand*{\figref}[1]{\figurename~\ref{#1}} \newcommand*{\figref}[1]{\figurename~\ref{#1}}
\captionsetup{labelfont={bf}} \captionsetup{labelfont={bf}}
@ -47,25 +67,25 @@
} }
} }
\abstract { \abstract {Physics engines are a fun and interesting way to learn
Physics engines are a fun and interesting way to learn about a lot of about a the laws of physics, as well as computer science. They
different subjects. First the theoretical concepts, such as the equations provide a real-time simulation of common physical phenomena, and
that dictate the motion of the objects, together with their components, need therefore illustrate theoretical concepts such as the equations that
to be thoroughly understood. Then there is the necessity of finding a way to dictate the motion of objects, The goal of this project was to
represent all of those concepts in a given programming language and to extend an existing physics engine built for demonstration purposes.
make them as efficient as possible so that the simulation runs fluidly. This engine was initially designed and developed to simulate
The task to be completed here was to extend an already existing circular objects (``balls'') in 2D. With this project, we intended
physics engine that only made circles bounce off each other. The extension to extend this engine to also simulate arbitrary polygons, again in
was focused on having the ability to generate some arbitrary polygons and a physically accurate way. The main technical challenges of the
make them bounce off each other in a physically accurate way. The main project is therefore the correct simulation of the dynamics of
issues that rose up during the development of the extension: determining the rigid, polygonal objects. In particular, we developed a model of
inertia of a arbitrary polygon, which is important for realistic polygonal rigid objects; we implemented a simulation of their
impacts; having an accurate collision detection system, which allows the inertial motion, possibly in the presence of a constant force field
engine to know when to make two polygons bounce off each other. Once those such as gravity; we detect collisions between objects; we compute
aspects were worked on and polished, the rest of the implementation went and then simulate the dynamic effects of collisions. The
smoothly. simulations are animated and displayed in real-time. It is also
therefore crucial that the simulation code be efficient to obtain
} smooth animations.}
\counterwithin{figure}{section} \counterwithin{figure}{section}
\counterwithin{equation}{section} \counterwithin{equation}{section}

BIN
figures/after.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 27 KiB

BIN
figures/before.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 38 KiB

BIN
figures/bounding_box.pdf Normal file

Binary file not shown.

View File

@ -3,7 +3,7 @@
\section{Calculations} \section{Calculations}
\label{appendix:calculations} \label{appendix:calculations}
\paragraph{Moment of inertia of rectangle} \subsection{Moment of inertia of rectangle}
\begin{equation} \begin{equation}
\label{eq:rect_moment_long} \label{eq:rect_moment_long}
\begin{split} \begin{split}
@ -22,7 +22,7 @@
\end{equation} \end{equation}
\newpage \newpage
\paragraph{Moment of inertia of sub-triangle of regular polygon} \subsection{Moment of inertia of sub-triangle of regular polygon}
Before starting the calculations, it is to be noted that according to Figure Before starting the calculations, it is to be noted that according to Figure
\ref{fig:subtriangle}, we have that \ref{fig:subtriangle}, we have that
$$ \tan\left(\frac{\theta}{2}\right) = \frac{\frac{l}{2}}{h} = \frac{l}{2h} $$ $$ \tan\left(\frac{\theta}{2}\right) = \frac{\frac{l}{2}}{h} = \frac{l}{2h} $$
@ -46,7 +46,7 @@ it will be useful to simplify the result of the integral.
\end{equation} \end{equation}
\newpage \newpage
\paragraph{Moment of inertia of sub-triangle of arbitrary polygon} Recall \subsection{Moment of inertia of sub-triangle of arbitrary polygon} Recall
equation \ref{eq:r} defines equation \ref{eq:r} defines
$$ \vec r = \alpha \vv{CA} + \beta \alpha \vv{AB} $$ $$ \vec r = \alpha \vv{CA} + \beta \alpha \vv{AB} $$
\begin{equation} \begin{equation}
@ -62,3 +62,37 @@ $$ \vec r = \alpha \vv{CA} + \beta \alpha \vv{AB} $$
&= \frac{\rho hb}{4} \left(\frac{1}{3}\vv{AB}^2 + \vv{AB} \cdot \vv{CA} + \vv{CA}^2\right) \\ &= \frac{\rho hb}{4} \left(\frac{1}{3}\vv{AB}^2 + \vv{AB} \cdot \vv{CA} + \vv{CA}^2\right) \\
\end{split} \end{split}
\end{equation} \end{equation}
\newpage
\subsection{Solving for impulse parameter}
\label{app:impulse_long}
We start with equation \ref{eq:vp2n}:
\begin{equation*}
\begin{split}
\vec v_{p2} \cdot \vec n &= - e \vec v_{p1} \cdot \vec n\\
\left(\vec v_{ap2} - \vec v_{bp2}\right)\cdot \vec n &= - e \vec v_{p1} \cdot \vec n\\
\left( \vec v_{a2} + \omega_{a2} \times \vec r_{ap} - \vec v_{b2} - \omega_{b2} \times \vec r_{bp}
\right)\cdot \vec n &= - e \vec v_{p1} \cdot \vec n\\
\end{split}
\end{equation*}
We now expand the bracket on the left-hand side using the equations \ref{eq:va2}
- \ref{eq:omega_b2} and then simplify with equation \ref{eq:vp1}.
\[
\begin{split}
\left( \vec v_{a1} + \frac{j\vec n}{m_a} + \omega_{a1} + \frac{\vec r_{ap} \times j\vec n}{I_a} \times \vec r_{ap} - \vec v_{b1} + \frac{j\vec n}{m_b} - \omega_{b1} + \frac{\vec r_{bp} \times j\vec n}{I_b} \times \vec r_{bp} \right)\cdot \vec n &= - e \vec v_{p1} \cdot \vec n\\
\left( \frac{j\vec n}{m_a} + \frac{\vec r_{ap} \times j\vec n}{I_a} \times \vec r_{ap} + \frac{j\vec n}{m_b} + \frac{\vec r_{bp} \times j\vec n}{I_b} \times \vec r_{bp} \right)\cdot \vec n &= - (1+e) \vec v_{p1} \cdot \vec n\\
\end{split}
\]
Using the triple scalar product rule, we can derive that
\[
j \left( \frac{1}{m_a} + \frac{1}{m_b} + \frac{\left(\vec r_{ap} \times \vec n\right)^2}{I_a} + \frac{\left( \vec r_{bp} \times \vec n \right)^2}{I_b}
\right) = -(1+e) \vec v_{p1} \cdot \vec n
\]
and therefore
\begin{equation}
j = \frac{ - (1+e) \cdot \vec v_{ap1} \cdot \vec n }{\frac{1}{m_a} + \frac{1}{m_b} +
\frac{\left( \vec r_{ap} \times \vec n \right)^2}{I_a} + \frac{\left( \vec
r_{bp} \times \vec n \right)^2}{I_b}}
\end{equation}

View File

@ -1 +1,26 @@
\section{Conclusion} \section{Conclusion}
The overall development of this extension was very instructive. Getting to
discover some theoretical concepts, having to understand them and then being
able to apply them in the code is a very rewarding experience. The objective of
project is considered accomplished. In figure \ref{fig:before_after}, we can see
a comparison of the state of the application before and after the extension.
\begin{figure}[H]
\centering
\hfill
\begin{subfigure}[]{.4\textwidth}
\centering
\includegraphics[width=.8\textwidth]{before}
\caption{Original state of the simulation}
\end{subfigure}
\hfill
\begin{subfigure}[]{.4\textwidth}
\centering
\includegraphics[width=.8\textwidth]{after}
\caption{State after the extension}
\end{subfigure}
\hfill\null
\caption{Before v. After}
\label{fig:before_after}
\end{figure}

View File

@ -1 +1,178 @@
\section{Implementation} \section{Implementation}
This section will dive into the actual implementing of the project, trying to
describe how the theoretical concepts in section \ref{sec:theory}.
\subsection{Structure}
The added code to the project can be divided into three parts
\begin{itemize}
\item polygon generator;
\item collision detection;
\item collision resolution.
\end{itemize}
Each part will be explained in the following sub-sections.
\subsubsection{Polygon generator}
Before talking about how we generate polygons, let's first talk about how the
polygons are defined in this project.
In the file \texttt{polygons.h} we define a polygon as
\begin{lstlisting}[caption={Polygon class (simplified)},label={lst:polygon}]
class polygon {
std::vector<vec2d> points;
vec2d center;
double angle;
double inertia;
double mass;
std::vector<vec2d> global_points();
vec2d speed;
double angular_speed;
}
\end{lstlisting}
Polygons possesses multiple fields. Firstly we have
\lstinline{std::vector<vec2d> points}, a collection of
\lstinline{vec2d} objects, which represent the set of ordered points that
compose the polygon. Those points are expressed in local coordinates, and the
center of mass of the polygon is placed the origin.
Now that we know how the polygon is composed, we want to move it around the
space it lives in, that is the purpose of the \lstinline{vec2d center} field.
It represents where the center of mass of the polygon is located in simulation.
Since the shapes also rotate, we need to keep track of the rotation of the
polygon, hence the use of \lstinline{double angle}. All of these fields are used
when computing the result of the method \lstinline{std::vector<vec2d> global_points()}.
Finally, we have \lstinline{double inertia} and \lstinline{double mass} which,
as their name suggest, store the value for the polygon's inertia and mass, that
are used to calculate the polygons final speed and angular speed, who are
represented by \lstinline{vec2d speed} and \lstinline{double angular_speed}.
\paragraph{Polygon generator} Now that we know how polygons are represented in
our simulation, we can generate them. In \texttt{polygon\_generator.h} (and its
implementation file \texttt{polygon\_generator.cc}), there are some functions
for that. In Listing \ref{lst:polygon_gen} we can see the signature of the
functions that generate
\begin{itemize}
\item a rectangle of a certain width and height;
\item a square of a certain side length;
\item a triangle given two side lengths and an angle between them;
\item a regular polygon;
\item an arbitrary polygon.
\end{itemize}
\begin{lstlisting}[caption={Polygon Generator header
file},label={lst:polygon_gen}]
namespace poly_generate {
polygon rectangle(double width, double height);
inline polygon square(double width) {
assert(width > 0);
return rectangle(width, width, label);
};
polygon triangle(double side1, double side2, double angle);
polygon regular(double radius, uint n_sides);
polygon general(std::vector<vec2d> points);
};
\end{lstlisting}
The implementation of those functions are fairly straight forward, they generate
the appropriate number of points in the correct place. After what they calculate
the mass and subsequent inertia of the polygon as shown in section
\ref{sub:moment}.
\subsubsection{Collision}
\label{sub:implementation-collision-detection}
The algorithm described in section~\ref{sub:vertex-collision} is implemented in
\texttt{collision.cc} and exposed to other modules with \texttt{collision.h}.
The module exposes a collision structure and a collides function.
\begin{lstlisting}[caption={Collision header file},label={lst:collision}]
struct collision {
bool collides = false;
vec2d impact_point;
vec2d n;
vec2d overlap;
};
extern collision collides(polygon& p, polygon& q);
\end{lstlisting}
The \lstinline{collides} function takes in two polygons and checks with
vertex-collision algorithm whether they collide or not. The result is return
through an instance of \lstinline{struct collision}. The \lstinline{bool collides}
and \lstinline{vec2d impact_point} fields is self-explanatory. The
\lstinline{vec2d n} is the normal vector that pushes \lstinline{polygon& p} away
from \lstinline{polygon& q}. Finally, \lstinline{vec2d overlap} is a scalar
multiplication of \lstinline{vec2d n}. Where the latter is a normalized vector,
the former represents how deep \lstinline{p} is in \lstinline{q}. If we push
\lstinline{p} the exact amount of \lstinline{overlap}, then \lstinline{p} would
just be touching \lstinline{q} with point \lstinline{impact_point}, and no
further overlap would occur.
In reality, we do not simply push \lstinline{p} by \lstinline{overlap}, but we
push \lstinline{p} and \lstinline{q} away from each other proportionally to
their mass.
\paragraph{Collision resolution} In \texttt{polygons.cc}, among other functions,
we apply the collision resolution, which simply uses the physics results found
in \ref{sub:resolution}.
\subsection{Optimization}
In order to optimise the collision detection and to avoid having to do a check
whether each vertex of every polygon was colliding with each edge of every other
polygon we decided to apply the bounding box acceleration structure.
The bounding box of an object is, in two dimensions, the rectangle that contains
the object and which the sides are aligned with axis of the coordinate system.
\begin{figure}[H]
\centering
\inputtikz[.7]{bounding_box}
\caption{Bounding box example}
\label{fig:bounding-box}
\end{figure}
To get the four coordinates that compose, we just perform a linear scan through
all the points that compose a polygon, record the minimum and maximum of both
they $x$ and $y$ coordinate.
Once those points have been found, we just perform some basic axis-aligned
rectangle collision detection, which is much less computationally expensive than
checking each vertex-edge pair.
The complexity of the collision detection algorithm went from $\mathcal{O}(n^2)$
to $\mathcal{O}(n)$, where $n$ is the total number of vertices across all
polygons in the simulation.
\subsection{Known issues}
The simulation has one major flaw, and it is inherent to the last part of
section \ref{sub:implementation-collision-detection} where we talk about the
\lstinline{overlap} vector: we said that this vector allows to displace both
polygons so that and the end of the collision resolution, they are not
overlapping anymore. This was done to avoid issues where the collision
resolution had taken place in one frame, but the polygons were still overlapping
in the following frame, which lead to polygons getting stuck together.
The unfortunate consequence of such a decision is that, when we start playing
around with the restitution coefficient and all the shapes start falling to the
ground, the shapes are not able to rest still. Since there is gravity, at each
frame their speed get updated to simulate its effect, but since they are lying
on the floor, the collision detection algorithm kicks in right away, resulting
in the polygon getting pushed completely up, since the ground has infinite mass.
This results in the polygon bouncing on the floor instead of lying down.

View File

@ -1,44 +1,61 @@
\section{Theoretical Background} \section{Theoretical Background}
\label{sec:theory}
The theoretical background is everything related to the physics part of the The theoretical background is everything related to the physics part of the
project. It covers the calculating the inertia of different types of polygons; project. It covers the calculation of the inertia of different types of polygons;
different algorithms to detect whether there is a collision between two different algorithms to detect whether there is a collision between
polygons; the resolution of the collision, i.e. finding the final any two polygons; and the resolution of the collision, that is, finding the final
velocity vectors and angular speed of those polygons. velocity vectors and angular speed of those polygons.
\subsection{Moment of inertia} \subsection{Moment of Inertia}
\label{sub:moment}
The inertia of an object refers to the tendency of an object to resist a change The inertia of an object refers to the tendency of an object to resist
of its state of motion or rest, it describes how the object behaves when forces a change of its state of motion or rest, it describes how the object
are applied to it. An object with a lot of inertia requires more force to change behaves when forces are applied to it. In particular, the mass of a
its motion, either to make it move if it's at rest or to stop it if it's already point object (inertial mass) determines the inertia of that objects.
moving. On the other hand, an object with less inertia is easier to set in Thus an object with a large mass requires more force to change its
motion or bring to a halt. motion, either to make it move from an initial state at rest or to
stop it if it is already moving. On the other hand, an object with a
smaller mass is easier to set in motion or to bring to a halt.
The moment of inertia is similar but is used in a slightly different context, it For objects that are not point-like, the motion includes a rotational
specifically refers to the rotational inertia of an object. It measures an component.
object's resistance to changes in its rotational motion and how its mass is
distributed with respect to is axis of rotation.
In the case of this project the axis of rotation is the one along the $z$-axis The moment of inertia is analogous to the mass in determining the
(perpendicular to the plane of the simulation) and placed at the barycenter of inertia of an object with respect to its rotational motion. In other
the polygon. words, the moment of inertia determines the rotational inertia of an
object, that is, the object's resistance to changes in its rotational
motion. The moment of inertia of an object is therefore computed with
respect to a rotation axis, and it is determined by the way the mass
of the object is distributed with respect to the chosen axis of
rotation.
The general formula for the moment of inertia is
For the purpose of this project, we are interested in a 2D simulation.
This means that the axis of rotation is always parallel to the
$z$-axis, perpendicular to the plane of the simulation. Furthermore,
we model the motion of an object using the barycenter of an object as
its reference point. We therefore compute the moment of inertia of a
polygon with respect to a rotation axis that is perpendicular to the
2D plane and that goes through the center of mass of the polygon.
In general, the moment of inertia about a rotation axis $Z$ of a point
mass $m$ at distance $r$ from $Z$ is $I_m=mr^2$. Therefore, for a
generic 2D object $Q$, we integrate over the whole 2D object:
\begin{equation} \begin{equation}
\label{eq:moment_general} \label{eq:moment_general}
I_Q = \int \vec r^2 \rho(\vec r) \diff \mathcal{A} I_Q = \int_{\mathcal{A}}\|\vec r\|^2 \rho(\vec r) \diff \mathcal{A}
\end{equation} \end{equation}
where $\rho$ is the density of object $Q$ in the point $\vec r$ across the where $\rho$ is the density of object $Q$ (mass per unit area) in the
small pieces of area $\mathcal A$ of the object. point $\vec r$ from the across the small pieces of area $\mathcal A$
of the object.
In our case, since we are implementing a 2D engine we can use the $\mathbb{R}^2$ In our case, we use 2D Cartesian coordinates centered in the point
coordinate systems, thus the formula becomes (axis) about which we compute he moment of inertia. Furthermore, we
$$ I_Q = \iint \rho(x, y) \vec r^2 \diff x\diff y$$ consider polygons of uniform density. We can therefore compute:
and since the requirements express that the mass of the polygons is spread
uniformly across its surface, the formula finally becomes
\begin{equation} \begin{equation}
\label{eq:moment} \label{eq:moment}
I_Q = \rho \iint x^2 + y^2 \diff x\diff y I_Q = \rho \iint(x^2 + y^2)\diff x\diff y
\end{equation} \end{equation}
The bounds of the integral depend on the shape of the polygon. In the following The bounds of the integral depend on the shape of the polygon. In the following
@ -264,6 +281,7 @@ to the overall polygon.
where, $P_{n+1} = P_1$ in the case of $i = n$. where, $P_{n+1} = P_1$ in the case of $i = n$.
\subsection{Collision detection} \subsection{Collision detection}
\label{sub:collision-detection}
Collision detection, as the name suggests, are the algorithms used to detect Collision detection, as the name suggests, are the algorithms used to detect
whether two polygons are colliding. The result of this procedure must be an whether two polygons are colliding. The result of this procedure must be an
@ -324,6 +342,7 @@ algorithm of our own. Moreover, SAT only supports convex polygons, which limits
the original objective of the project, which was to have any arbitrary polygon. the original objective of the project, which was to have any arbitrary polygon.
\subsubsection{Vertex collisions} \subsubsection{Vertex collisions}
\label{sub:vertex-collision}
The solution that was adopted for the project, after trying SAT, was a more The solution that was adopted for the project, after trying SAT, was a more
intuitive one, developed by Prof. Carzaniga. The idea is simple: check if a intuitive one, developed by Prof. Carzaniga. The idea is simple: check if a
@ -462,3 +481,164 @@ find the normal. The results look realistic enough to be accepted.
\subsection{Collision resolution} \subsection{Collision resolution}
\label{sub:resolution} \label{sub:resolution}
The collision resolution is the last step in the processing of the collision.
Algorithmically, it is much less heavy than collision detection, since, once the
simulation has two colliding polygons, a point of impact and a normal vector,
it's just a case of applying the rigid body physics formulas to the polygons
that are colliding. This part has been helped a lot by the works of Erik Neumann
\cite{collision:resolution-site} and Chris Hecker
\cite{collision:resolution-paper}.
\begin{figure}[H]
\centering
\inputtikz[.5]{collision_resolution}
\caption{Collision resolution between polygons $A$ and $B$}
\label{fig:collision_resolution}
\end{figure}
\paragraph{Variable definition} Before getting into any maths, let's define some
variables that we are going to use
\begin{multicols}{2}
\begin{itemize}
\item $m_a, m_b = $ mass of the bodies $A$ and $B$
\item $\vec r_{ap} = $ distance vector from center of mass of body $A$ to
point $P$
\item $\vec r_{bp} = $ distance vector from center of mass of body $A$ to
point $P$
\item $\omega_{a1}, \omega_{b1} = $ initial angular velocity of bodies
$A, B$
\item $\omega_{a2}, \omega_{b2} = $ final angular velocity of bodies
$A, B$
\item $\vec v_{a1}, \vec v_{b1} =$ initial velocities of center of mass
bodies $A, B$
\item $\vec v_{a2}, \vec v_{b2} =$ final velocities of center of mass
bodies $A, B$
\item $\vec v_{ap1}=$ initial velocity of impact point $P$ on body $A$
\item $\vec v_{bp1}=$ initial velocity of impact point $P$ on body $B$
\item $\vec v_{p1}=$ initial relative velocity of impact points on body $A, B$
\item $\vec v_{p2}=$ final relative velocity of impact points on body $A, B$
\item $\vec n=$ normal vector
\item $ e=$ elastic coefficient (0 = inelastic, 1 = perfectly elastic)
\end{itemize}
\end{multicols}
Some of those variables, like the ones in the left column, are already given by
the simulation, some of them have been computed thanks to the algorithms in
section \ref{sub:collision-detection} (normal vector and position of point $P$),
and some have still to be defined mathematically, such as $\vec v_{ap1},\vec
v_{bp1},\vec v_{p1}$ and $\vec v_{p2}$. So the velocity vectors of impact
point $P$ on both bodies, before the collision, are
\begin{equation} \label{eq:vabp1}
\begin{split}
\vec v_{ap1} = \vec v_{a1} + \omega_{a1} \times \vec r_{ap} \\
\vec v_{bp1} = \vec v_{b1} + \omega_{b1} \times \vec r_{bp}
\end{split}
\end{equation}
Similarly, the final velocities are
\begin{equation} \label{eq:vabp2}
\begin{split}
\vec v_{ap2} = \vec v_{a2} + \omega_{a2} \times \vec r_{ap}\\
\vec v_{bp2} = \vec v_{b2} + \omega_{b2} \times \vec r_{bp}
\end{split}
\end{equation}
Here we are regarding the angular velocity as a 3-dimensional vector
perpendicular to the plane, so that the cross product is calculated as
$$ \omega \times \vec r = \begin{pmatrix} 0\\0\\\omega \end{pmatrix} \times
\begin{pmatrix} r_x\\r_y\\0 \end{pmatrix} = \begin{pmatrix} -\omega r_y \\
\omega r_x \\0\end{pmatrix} $$
We these variables, we can finally define the relative velocities $\vec
v_{p1}$ and $\vec v_{p2}$
\begin{equation}
\label{eq:vp1}
\begin{split}
\vec v_{p1} = \vec v_{ap1} - \vec v_{bp2}\\
\vec v_{p2} = \vec v_{ap2} - \vec v_{bp2}
\end{split}
\end{equation}
If we expand by using \ref{eq:vabp1} and \ref{eq:vabp2}, we get
\begin{equation}
\begin{split}
\vec v_{p1} = \vec v_{a1} + \omega_{a1} \times \vec r_{ap} - \vec v_{b1}
+ \omega_{b1} \times \vec r_{bp}\\
\vec v_{p2} = \vec v_{a2} + \omega_{a2} \times \vec r_{ap} - \vec v_{b2} + \omega_{b2} \times \vec r_{bp}
\end{split}
\end{equation}
The relative velocity of point $P$ along the normal vector $\vec n$ is
$$ \vec v_{p1} \cdot \vec n $$
Note that for a collision to occur this relative normal velocity must be
negative (that is, the objects must be approaching each other). Let e be the
elasticity of the collision, having a value between 0 (inelastic) and 1
(perfectly elastic). We now make an important assumption in the form of the
following relation
\begin{equation}
\label{eq:vp2n}
\vec v_{p2} \cdot \vec n = - e \vec v_{p1} \cdot \vec n
\end{equation}
This says that the velocity at which the objects fly apart is proportional to
the velocity with which they were coming together. The proportionality factor is
the elasticity $e$.
\paragraph{Collision Impulse} In simple terms, collision impulse refers to the
change in momentum experienced by an object during a collision. It is a measure
of the force applied to an object over a short period of time. We imagine that
during the collision there is a very large force acting for a very brief period
of time. If you integrate (sum) that force over that brief time, you get the
impulse.
In our simulation, we are assuming that no friction is happening during the
collision, so that the only force we have to consider is the one of along the
normal vector $\vec n$. The friction would create a force perpendicular to the
normal, and it would make things a little too complicated for the scope of this
project.
Since the only force we consider is the normal one, we can consider the net
impulse to be $j \vec n$, where $j$ is the impulse parameter. The body $A$ will
experience the net impulse $j \vec n$ and body $B$ will experience it's negation
$- j \vec n$ since the force that $B$ experience is equal an opposite to the one
experienced by $A$. The impulse is a change in momentum. Momentum has units of
velocity times mass, so if we divide the impulse by the mass we get the change
in velocity. We can relate initial and final velocities as
\begin{equation}
\label{eq:va2}
\vec v_{a2} = \vec v_{a1} + \frac{j\vec n}{m_a}
\end{equation}
\begin{equation}
\label{eq:vb2}
\vec v_{b2} = \vec v_{b1} - \frac{j\vec n}{m_b}
\end{equation}
For the final angular speed, it's change from the impulse $j\vec n$ is given by
$\vec r_{ap} \times j \vec n$. We can divide the result by the moment of
inertia, which was calculated in section \ref{sub:moment}, to get convert the
change in angular momentum into change in angular speed. Similarly as above, we
can relate the initial and final angular velocity as
\begin{equation}
\label{eq:omega_a2}
\omega_{a2} = \omega_{a1} + \frac{\vec r_{ap} \times j\vec n}{I_a}
\end{equation}
\begin{equation}
\label{eq:omega_b2}
\omega_{b2} = \omega_{b1} - \frac{\vec r_{bp} \times j\vec n}{I_b}
\end{equation}
\paragraph{Solving for the impulse parameter}
Now we have everything we need to solve for the impulse parameter $j$. The bulk
of the calculations can be found in section \ref{app:impulse_long} of the
appendix \ref{appendix:calculations}. But basically we start with equation
\ref{eq:vp2n} and we end up with
\begin{equation}
\label{eq:j}
j = \frac{ - (1+e) \cdot \vec v_{ap1} \cdot \vec n }{\frac{1}{m_a} + \frac{1}{m_b} +
\frac{\left( \vec r_{ap} \times \vec n \right)^2}{I_a} + \frac{\left( \vec
r_{bp} \times \vec n \right)^2}{I_b}}
\end{equation}

View File

@ -37,7 +37,7 @@
\draw [style=none] (15.center) to (14.center); \draw [style=none] (15.center) to (14.center);
\draw [style=none] (14.center) to (13.center); \draw [style=none] (14.center) to (13.center);
\draw [style=none] (13.center) to (12.center); \draw [style=none] (13.center) to (12.center);
\draw [style=Axis] (0.center) to (16); \draw [style=Vector] (0.center) to (16);
\draw [style=Dotted] (0.center) to (21); \draw [style=Dotted] (0.center) to (21);
\draw (23.center) to (26.center); \draw (23.center) to (26.center);
\draw (26.center) to (25.center); \draw (26.center) to (25.center);

71
tikzs/bounding_box.tikz Normal file
View File

@ -0,0 +1,71 @@
\begin{tikzpicture}
\begin{pgfonlayer}{nodelayer}
\node [style=none] (0) at (0, -1) {};
\node [style=none] (1) at (0, 8) {};
\node [style=none] (2) at (-1, 0) {};
\node [style=none] (3) at (15, 0) {};
\node [style=none] (4) at (4, 4) {};
\node [style=none] (5) at (2, 1) {};
\node [style=none] (6) at (6, 3) {};
\node [style=none] (7) at (9, 3) {};
\node [style=none] (8) at (11, 2) {};
\node [style=none] (9) at (10.5, 6.5) {};
\node [style=none] (10) at (12.75, 5.75) {};
\node [style=none] (11) at (13.5, 3.5) {};
\node [style=none] (12) at (8, 5) {};
\node [style=none] (13) at (2, 4) {};
\node [style=none] (14) at (6, 1) {};
\node [style=none] (15) at (6, 4) {};
\node [style=none] (16) at (8, 2) {};
\node [style=none] (17) at (8, 6.5) {};
\node [style=none] (18) at (13.5, 6.5) {};
\node [style=none] (19) at (13.5, 2) {};
\node [style=none] (20) at (2, 0) {};
\node [style=none] (21) at (6, 0) {};
\node [style=none] (22) at (0, 4) {};
\node [style=none] (23) at (0, 1) {};
\node [style=none] (24) at (-0.25, 1) {$A_{y_{\min}}$};
\node [style=none] (25) at (-0.25, 4) {$A_{y_{\max}}$};
\node [style=none] (26) at (2, -0.25) {$A_{x_{\min}}$};
\node [style=none] (27) at (6, -0.25) {$A_{x_{\max}}$};
\node [style=none] (28) at (8, 0) {};
\node [style=none] (29) at (8, -0.25) {$B_{x_{\min}}$};
\node [style=none] (30) at (13.5, 0) {};
\node [style=none] (31) at (13.5, -0.25) {$B_{x_{\max}}$};
\node [style=none] (32) at (0, 2) {};
\node [style=none] (33) at (-0.25, 2) {$B_{y_{\min}}$};
\node [style=none] (34) at (0, 6.5) {};
\node [style=none] (35) at (-0.25, 6.5) {$B_{y_{\max}}$};
\node [style=none] (36) at (4, 2.75) {$A$};
\node [style=none] (37) at (10.5, 4.5) {$B$};
\end{pgfonlayer}
\begin{pgfonlayer}{edgelayer}
\draw [style=Vector] (0.center) to (1.center);
\draw [style=Vector] (2.center) to (3.center);
\draw (5.center) to (4.center);
\draw (4.center) to (6.center);
\draw (6.center) to (5.center);
\draw (12.center) to (9.center);
\draw (9.center) to (10.center);
\draw (10.center) to (11.center);
\draw (11.center) to (8.center);
\draw (8.center) to (7.center);
\draw (7.center) to (12.center);
\draw [style=Dashed] (13.center) to (5.center);
\draw [style=Dashed] (5.center) to (14.center);
\draw [style=Dashed] (14.center) to (15.center);
\draw [style=Dashed] (15.center) to (13.center);
\draw [style=Dashed] (16.center) to (17.center);
\draw [style=Dashed] (17.center) to (18.center);
\draw [style=Dashed] (18.center) to (19.center);
\draw [style=Dashed] (19.center) to (16.center);
\draw [style=Dotted] (19.center) to (30.center);
\draw [style=Dotted] (16.center) to (28.center);
\draw [style=Dotted] (16.center) to (32.center);
\draw [style=Dotted] (13.center) to (22.center);
\draw [style=Dotted] (17.center) to (34.center);
\draw [style=Dotted] (5.center) to (20.center);
\draw [style=Dotted] (5.center) to (23.center);
\draw [style=Dotted] (14.center) to (21.center);
\end{pgfonlayer}
\end{tikzpicture}

View File

@ -17,8 +17,8 @@
\node [style=none] (14) at (2, .5) {\Large$\frac{\theta}{2}$}; \node [style=none] (14) at (2, .5) {\Large$\frac{\theta}{2}$};
\end{pgfonlayer} \end{pgfonlayer}
\begin{pgfonlayer}{edgelayer} \begin{pgfonlayer}{edgelayer}
\draw [style=Axis] (0.center) to (1.center); \draw [style=Vector] (0.center) to (1.center);
\draw [style=Axis] (2.center) to (3.center); \draw [style=Vector] (2.center) to (3.center);
\draw (4.center) to (5.center); \draw (4.center) to (5.center);
\draw (5.center) to (6.center); \draw (5.center) to (6.center);
\draw (6.center) to (4.center); \draw (6.center) to (4.center);

View File

@ -25,7 +25,7 @@
\draw (1.center) to (2.center); \draw (1.center) to (2.center);
\draw (2.center) to (3.center); \draw (2.center) to (3.center);
\draw (3.center) to (0.center); \draw (3.center) to (0.center);
\draw [style=Axis] (16.center) to (11.center); \draw [style=Vector] (16.center) to (11.center);
\draw [style=Axis] (15.center) to (12.center); \draw [style=Vector] (15.center) to (12.center);
\end{pgfonlayer} \end{pgfonlayer}
\end{tikzpicture} \end{tikzpicture}

View File

@ -25,8 +25,8 @@
\draw (1.center) to (2.center); \draw (1.center) to (2.center);
\draw (2.center) to (3.center); \draw (2.center) to (3.center);
\draw (3.center) to (0.center); \draw (3.center) to (0.center);
\draw [style=Axis] (16.center) to (11.center); \draw [style=Vector] (16.center) to (11.center);
\draw [style=Axis] (15.center) to (12.center); \draw [style=Vector] (15.center) to (12.center);
\draw [style=Axis] (21.center) to (22.center); \draw [style=Vector] (21.center) to (22.center);
\end{pgfonlayer} \end{pgfonlayer}
\end{tikzpicture} \end{tikzpicture}