bAST 2020 VerkadeT

Rijksuniversiteit Groningen
Bachelor Thesis
Simulating Galaxy Collisions in Python

for Astronomy Education
Author: Supervisor:
Thijs Verkade (s3479838) Dr. Jake Noel-Storr
July 3, 2020
Thijs Verkade
Abstract
In this thesis, the process behind producing an application in python which simulates galaxy
collisions is explained. The aim of the application is for it to be used by students in order to
introduce them to the concepts of galaxy collisions in an educational manner. The application
is based upon a java application called GalCrash, created by Chris Mihos. The application is
tested by students and the educational value of the application is then discussed based on the
feedback provided.
i
Thijs Verkade
Acknowledgements
I would like to thank Dr. Jake Noel-Storr for supplying the thesis topic, as well as for their support and
guidance in producing the application, as well as organizing the evaluation of the application. I would
also like to thank the group of students, along with the teaching assistant, who tested the application,
and provided valuable feedback that was used to improve the application.
ii
Thijs Verkade CONTENTS
Contents
1 Introduction 1
2 Theory 2
3 Coding 4
3.1 Object-Oriented Programming . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
3.1.1 Classes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
3.1.2 Instances . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
3.2 Implementing OOP in Python . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
3.2.1 Galaxy Class . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
3.2.2 StarGalaxy Class . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
3.2.3 Orbit Class . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
3.3 Creating the Graphical User Interface . . . . . . . . . . . . . . . . . . . . . . . . 14
3.4 Combining the code with the GUI . . . . . . . . . . . . . . . . . . . . . . . . . . 16
3.4.1 Initializing the galaxies and stars . . . . . . . . . . . . . . . . . . . . . . . 16
3.4.2 Evolving the galaxies and stars over time . . . . . . . . . . . . . . . . . . 22
3.4.3 Interrupting the simulation . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.4.4 Running the Application . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
4 Evaluating the educational value of the Application 27

4.1 Discussing the feedback from tutorial session . . . . . . . . . . . . . . . . . . . . 27
4.2 Changes made to the application . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
4.3 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
5 Conclusion 30
5.1 Further Improvements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
A Derivation of Rotation Matrix 31
B Evaluation Report 32
iii
Thijs Verkade 1. Introduction
1 Introduction
Studying the interactions between galaxies play a fundamental role in understanding how galaxies
have evolved into the shapes and sizes that are present today. The colliding and merging of
galaxies was a much more common process in the Early Universe, where astronomers found that
very young galaxies often resembled starburst galaxies with multiple nuclei and peculiar shapes
(Fraknoi et al., 2016). Merging galaxies tend to induce massive bursts of star formation and
completely reshape the galaxies on the scale of millions of years. To be able to understand the
properties of current galaxies and how these came about, it is thus extremely important to study
these galaxy interactions to understand how these changes occur. As galaxy interactions occur
over a timescale of millions of years, through use of computer simulations, one can study the
whole process of a galaxy collision or a merger without having to live an eternity.
Roughly 20 years ago, such an application was created by Chris Mihos that could run these
types of simulations. The application, named GalCrash, can simulate the interaction of two
spiral galaxies. In this application, the user can change a few of the properties of the galaxies
such as their relative masses, their separation, and the orientation of each galaxy’s disk. By
altering these properties, one can create a variety of different simulations that produce different
results to show the impact that each of these parameters have on each respective galaxy and
the resultant merger. This can then also be used to recreate collisions between real interacting
galaxies that we can observe in the sky. The code for this application was written in Java, and
is fairly outdated, while also lacking information in terms of how the code is written. It gives
little information to the reader as to what processes are being carried out in the code and why
this is being done.
The goal of this project, is to re-write the GalCrash application in Python, and do this in such
a way that it can be used as an education tool to teach students about galaxy interactions. By
updating the GalCrash application to be run in Python, the code also becomes more accessible
to students studying Astronomy, as Python is a much more common language being taught to
students studying Astronomy. In this way, students can still simulate galaxy interactions, while
also gaining more insight into how such an application is created. By writing the code in such a
way that the steps carried out are described, it becomes much easier for students to understand
the purpose of each line of code. The graphical user interface (GUI) for the Python application
written was created using Qt Designer, a tool for designing and building GUIs with widgets.
In this report, we will discuss the physics behind the galaxy interactions being simulated, as
well as walk-through the coding behind the simulation step-by-step, explaining what processes
were used to develop the final application. The educational value of the application will also be
evaluated, looking at the accessibility of the application and whether or not students find it to
be a valuable resource. Further improvements and features that could be added in the future
will also be discussed.
1
Thijs Verkade 2. Theory
2 Theory
The simulations carried out in this project are done between two spiral galaxies. Spiral galaxies
contain a disk of stars and gas arranged in a spiral pattern. Spiral galaxies usually contain a
central bulge that resembles a small elliptical galaxy. The bulge consists of older stars while
the disk is usually made up of a large range of ages, partly young stars (Gallagher & Sparke,
2007). According to a 2010 Hubble Space Telescope survey, around 70% of the galaxies observed
were classified as spirals (Delgado-Serrano et al., 2010). Therefore, when looking at galaxy
interactions, it makes sense to study the interaction between two spiral galaxies, as these have
the highest likelihood to interact with one another due to the fact that they are so common.
The process of star formation within galaxies are linked to the growth of galaxies through galaxy
interactions. The most common type of interaction between galaxies are that of a galaxy with
a minor companion, due to the greater fractional abundance of low luminosity galaxies (Lambas
et al., 2012). This is whats known as a minor merger. In this scenario, when the galaxies merge,
the larger galaxy remains relatively unchanged, while the small galaxy is completely disrupted
or its core comes to rest at the center of the larger galaxy.
A major merger occurs when the merging galaxies have similar masses. In such a case, due
to the violently changing gravitational fields, the merger remnant does not resemble either of
the original galaxies very well. In a major merger, the galaxies merge into a single steady-state
system.
In both minor and major mergers, dynamical friction plays an important role in how galaxies
collide. Dynamical friction is a pseudo-frictional force which occurs when a massive object moves
through a sea of lower mass particles (Mihos et al., 1999). Due to the large gravitational force
that the massive object is exerting on the lower mass particles, a high concentration of smaller
particles start to form behind the massive object. This large concentration of small particles
behind the massive object then exerts a collective gravitational force on the massive object,
causing it to decelerate. This is illustrated in Figure 1 below.
(a) (b) (c)
Figure 1: (a), (b) and (c) show the effect of dynamical friction on a massive object due to the
increasing density of the sea of stars and dark matter forming behind the object. The gray
represents a low density region while the black represents a high density region. (“Dynamical
Friction”, n.d.)
2
Thijs Verkade 2. Theory
Dynamical friction causes the galaxies to brake on their orbit and then merge. In our simulation,
dynamical friction does not occur, due to the simplicity of our simulation. However, this friction
can be analytically approximated using the Chandrasekhar formula. The Chandrasekhar formula
states that galaxies should feel a frictional force due to dynamical friction proportional to the
mass of the object times the density of the background sea divided by the velocity of the object
squared (Mihos et al., 1999) (1).
Mρ
F ∝ (1)
v2
To understand this process, we look at the case of a minor merger. We consider a companion
(small) galaxy with mass Mc , travelling through a population of field stars of individual mass ma ,
where ma << Mc . In this derivation, we assume that the companion galaxy acts as a point mass,
even though we know that this is not the case. The field stars present are part of a much larger
host galaxy of mass Mg where Mg >> Mc . Here, the assumption is taken that the mass of the
host galaxy is so much larger that it can be approximated as infinite and homogeneous. In such
a scenario, then the dominant effect of the encounters is dynamical friction, which decelerates
the companion galaxy at a rate of (Aceves & Colosimo, 2006):
~vc − va
Z
dvc
= −4πG Mc ma ln Λ d3 va f (va )
2
(2)
dt |~vc − va |3
Then, if the field stars have an isotropic velocity distribution, then equation 2 yields a simpler
expression for the dynamical friction, namely Chandrasekhar’s dynamical friction formula:
Z vc
dvc 2 2 2 ~vc
= −16π G Mc ma ln Λ dva va f (va ) (3)
dt 0 vc3
In the case for which vc is sufficiently large, then our integral in 3 converges to a definite limit
equal to the number density n divided by 4π. This then provides us with the following formula
for dynamical friction:
dvc ~vc
= −4πG2 Mc ma na ln Λ 3 (4)
dt vc
This can then be written in terms of the density ρa as follows,
dvc ~vc
= −4πG2 Mc ρa ln Λ 3 (5)
dt vc
This exact formula is not used in the simulation, but a very similar formula is used which mimics
the Chandrasekhar expression is used. This expression can be found in the following section of
the report.
3
Thijs Verkade 3. Coding
3 Coding
3.1 Object-Oriented Programming
To generate the code of the new GalCrash application, a programming paradigm called Object-
Oriented Programming (OOP) was used. OOP is an approach primarily used for modeling
real-world things, where real-world entities are modeled as software objects. For this project,
we have several different types of objects, each with their own properties. An example of an
object present in this project is a star, which has the following properties: mass, position and
velocity. For each such an object, as well as having properties, they also have functions which
can be carried out, such as "MoveStar", which requires a time-step input, and will provide new
positions and velocities so we can see the evolution of the properties of the star through time.
3.1.1 Classes
To be able to achieve this, we create whats called a class in Python. A class is essentially a
blueprint for how something should be defined. Classes are generally used when you want to
represent something much more complicated, with many different properties. Such classes can
thus track properties about an object in an organized matter. Continuing with our example of
stars, we created a class called star, in which we would like to know the mass of the star, its
position in three dimensions, as well as its velocity in three dimensions. It is important to note
that classes do not provide any real content themselves, and so even though we have defined
specifications for a star, stating that a mass, position and velocity are necessary for defining a
star, it will not actually state these properties.
3.1.2 Instances
With these classes, we use python objects called Instances. Simply put, instances are copies of
the class with actual values. Instead of being a blueprint, an instance is now a fully realized
object belonging to a specific class. You can create multiple instances, each with different values,
but all containing the same information. This allows us to create multiple stars, each with their
own unique properties, in such a manner that these properties can easily be called upon and
used. Now that we understand the basics of OOP, and know what a class is and what an object
of a class is, we must then look at how this is implemented in Python.
3.2 Implementing OOP in Python

3.2.1 Galaxy Class
To look at how we implement classes, let us take a look at an example from our application.
In this particular example, we look at our Galaxy class, which includes several functions which
alter our parameters.
4
Thijs Verkade 3.2 Implementing OOP in Python
c l a s s Galaxy :
2 " " "A c l a s s used t o d e f i n e t h e i n i t i a l p a r a m e t e r s o f a g e n e r a t e d g a l a x y . I t
a l s o c o n t a i n s f u n c t i o n s t h a t a r e used i n c a l c u l a t i n g t h e e v o l u t i o n o f t h e
g a l a x y . Here we g e n e r a l l y d e a l with t h e c e n t e r o f t h e g a l a x y when l o o k i n g a t
p a r a m e t e r s such a s p o s i t i o n . " " "
d e f __init__ ( s e l f , galmass , ahalo , vhalo , r t h a l o , g a l p o s , g a l v e l ) :
4 s e l f . galmass = galmass
s e l f . ahalo = ahalo
6 s e l f . vhalo = vhalo
s e l f . rthalo = rthalo
8 s e l f . galpos = galpos
self . galvel = galvel
10 s e l f . g a l a c c = np . f u l l ( ( 3 , 1 ) , 0 . )
Listing 1: Creating a class with attributes
The first thing we do when creating a class is to create attributes. Using the init() method,
we specify an object’s initial attributes and give them their default state. Our method has
several attributes, as can be seen from Listing 1. galmass represents the galaxy mass, ahalo is
the acceleration of the dark matter halo, vhalo is the velocity of the dark matter halo, rthalo
is the radius of the dark matter halo, galpos and galvel are the three dimensional position and
velocity of the galaxy respectively. Along with these attributes, a self variable is also required.
The self variable refers to the object itself. The self variable is also an instance of our galaxy class,
and instances of a class have varying values, allowing us to assign different values to different
instances. The self variable helps keep track of individual instances of a class, allowing us to
specify different parameters for each individual instance, such that we can have two galaxies with
different properties.
Having now defined the basic blueprint for a galaxy, we can begin to perform operations with the
attributes of our objects. To do this, we use so called instance methods, which is very similar to
a function in Python. The first argument of an instance method is always self. A simple example
of an instance method can be found in Listing 2. Looking at our instance method setPosVel,
we can see that it simply changes the position and the velocity of the galaxy depending on the
input. This particular instance is only used when initializing our system, as it gives the inital
positions and velocities of the particular galaxy.
d e f s e t P o s v e l ( s e l f , pos , v e l ) :
2 s e l f . g a l p o s = pos
s e l f . galvel = vel
Listing 2: Creating Instance Methods to perform operations on our galaxy
The next instance method used is called scaleMass and can be found in Listing 3. This instance
is used to define the parameters of the companion galaxy, based on the input which is the ratio
of the mass of the companion galaxy to the main galaxy. Both galaxies are initially created with
the default parameters from Chris Mihos’ java version. So essentially, we start by creating two
identical galaxies with the same properties. Then, based on the input of the ratio of masses given
by the user, we scale the parameters of the companion galaxy based on the scaleMass instance
method. If the input of the mass ratio is given as 1 (default), then all parameters remain identical.
Otherwise, the mass, ahalo, vhalo and rthalo all change according to the formula provided in the
code.
5
2 d e f s c a l e M a s s ( s e l f , massFact ) :
" " "The S c a l e Mass Function c a l c u l a t e s t h e p a r a m e t e r s o f t h e companion
g a l a x y . This i s done by t a k i n g t h e r a t i o o f t h e mass o f t h e companion g a l a x y
t o t h e main g a l a x y a s t h e i n p u t . " " "
4 s e l f . g a l m a s s = s e l f . g a l m a s s ∗ massFact
s e l f . v h a l o = 1 . 0 ∗ massFact ∗ ∗ 0 . 2 5
6 s e l f . a h a l o = 0 . 1 ∗ massFact ∗ ∗ 0 . 5
a2 = − s e l f . g a l m a s s / ( s e l f . v h a l o ∗ ∗ 2 )
8 a1 = −2.∗ s e l f . a h a l o ∗ s e l f . g a l m a s s / ( s e l f . v h a l o ∗ ∗ 2 )
a0 = − s e l f . g a l m a s s ∗ ( s e l f . a h a l o ∗ ∗ 2 ) / ( s e l f . v h a l o ∗ ∗ 2 )
10 q = a1 / 3 . 0 − ( a2 ∗ ∗ 2 ) / 9 . 0
r = ( a1 ∗ a2 −3.∗ a0 ) / 6 . 0 − ( a2 ∗ ∗ 3 ) / 2 7 . 0
12
s 1 = ( r + np . s q r t ( ( q ∗ ∗ 3 ) +( r ∗ ∗ 2 ) ) ) ∗ ∗ 0 . 3 3 3
14 s 2 = ( ( r − np . s q r t ( ( q ∗ ∗ 3 ) +( r ∗ ∗ 2 ) ) ) ∗ ∗ 0 . 3 3 3 )
16 s e l f . r t h a l o =( s 1+s 2 )−a2 /3
Listing 3: Instance Method "Scale Mass" used to scale galaxy parameters based on user input
Another instance method used is namely the MoveGalaxy method shown in listing 4. As the
name suggests, the MoveGalaxy method simply advances our galaxy forward in a given timestep
dt. To be able to do this, the galaxies position, velocity and acceleration must be known. The
positions and velocities of the galaxies are well known as these are attributes of the galaxy. The
acceleration, however, needs to be calculated.
d e f MoveGalaxy ( s e l f , dtime ) :
2 " " "Move Galaxy e v o l v e s t h e p o s i t i o n and v e l o c i t y o f t h e g a l a x y with i n p u t
dtime " " "
newpos = s e l f . g a l p o s + s e l f . g a l v e l ∗ dtime + 0 . 5 ∗ s e l f . g a l a c c ∗ ( dtime ∗ ∗ 2 )
4 newvel = s e l f . g a l v e l + s e l f . g a l a c c ∗ dtime
6 s e l f . g a l p o s = newpos ;
s e l f . g a l v e l = newvel ;
Listing 4: Instance method "Move Galaxy" which evolves the parameters of the galaxy as time
progresses"
The acceleration can be calculated using Newton’s well-known law of universal gravitation,
in which we have that an object at a distance r from an object of mass M will fell an acceleration
due to gravity equivalent to:
−GM
g= (6)
r2
Using this method, we are treating our galaxies as if they have fixed shapes, when in reality, the
shape of the galaxy is constantly changing as it interacts gravitationally with the other galaxy.
If we wanted to calculate the acceleration in a self-consistent matter, then we would need to
calculate all the gravitational forces between all the stars and dark matter particles which make
up the galaxies, which is extremely expensive computationally, and therefore not feasible for a
simulation of this nature. Having said this however, although we sacrifice some accuracy, using
the Chandrasekhar dynamical friction formula (discussed in previous section), we can allow the
galaxies to merge in a way that more closely resembles reality. We will come back to this concept
later.
6
To calculate the acceleration, we have thus introduced an instance method with an input posin
(see listing 5. posin is the position of the interacting object. For example, if we want to calculate
the acceleration of our main galaxy M with respect to our companion galaxy C, then we would use
the position of C as an input. From here, our method calculates the difference in position between
our two galaxies, and uses this to calculate the gravitational acceleration from equation 6. In
this method, we use a different instance method to calculate the mass called InteriorMass. How
the interior mass is calculated within this instance method is discussed in the next paragraph.
The next instance method calculates the gravitational potential energy (Listing 6) in a very
similar manner to how the acceleration is calculated, the only difference being the formula for
the gravitational potential 7
GM
U =− (7)
r
2
def Acceleration ( s e l f , posin ) :
4 """ A c c e l e r a t i o n f u n c t i o n takes the p o s i t i o n o f the o b j e c t ( s t a r / galaxy ) as
i n p u t and r e t u r n s t h e a c c e l e r a t i o n o f t h e o b j e c t . " " "
G = 1.0
6 dpos = p o s i n − s e l f . g a l p o s
#dx = dpos [ 0 ]
8 #dy = dpos [ 1 ]
#dz = dpos [ 2 ]
10 r = np . s q r t ( np . sum ( dpos ∗ ∗ 2 , a x i s = 0 ) )
#r = np . s q r t ( ( dx ∗ ∗ 2 ) +(dy ∗ ∗ 2 ) +(dz ∗ ∗ 2 ) )
12 AccMag= −(G∗ s e l f . I n t e r i o r M a s s ( r ) ) / ( r ∗ ∗ 2 )
c a l c a c c = ( dpos ∗AccMag ) / r
14
return calcacc
Listing 5: Instance method "Acceleration" used to calculate acceleration of galaxy
2 def Potential ( s e l f , posin ) :

" " " P o t e n t i a l c a l c u l a t e s t h e p o t e n t i a l o f t h e o b j e c t ( s t a r / g a l a x y ) based on
the input p o s i t i o n given . """
4 G= 1 . 0 ;
6 dpos = p o s i n − s e l f . g a l p o s
#dx = dpos [ 0 ]
8 #dy = dpos [ 1 ]
#dz = dpos [ 2 ]
10 r = np . s q r t ( np . sum ( dpos ∗ ∗ 2 , a x i s = 0 ) )
#r = np . s q r t ( ( dx ∗ ∗ 2 ) +(dy ∗ ∗ 2 ) +(dz ∗ ∗ 2 ) )
12 pot= G∗ s e l f . I n t e r i o r M a s s ( r ) / r
14 r e t u r n pot
Listing 6: Instance method "Potential" used to calculate the potential energy of our galaxy
The next instance method has already been mentioned previously, namely the InteriorMass
instance method (Listing 7. This method returns the interior mass of the object dependent
on the input r. r is also used in the potential and acceleration methods and is defined as the
magnitude of the differential vector of the position of the star/galaxy with respect to the galaxy
being called upon. From here, we simply have two formulas for calculating the interior mass,
dependent on whether or not the magnitude r lies within the dark matter halo radius of the
7
galaxy or not. If we consider the main galaxy M and the companion galaxy C again, if we wish
to calculate the interior mass of our main galaxy, then r is the magnitude of the difference in
position between C and M. Then, if this magnitude is less than the dark matter halo radius of
M, this means that our companion galaxy lies within the dark matter halo radius of M. If this
is the case, then the interior mass of our galaxy M is calculated using equation 8.
2 × r3
vH
MI = (8)
(aH + r)2
Otherwise, if we have that r is greater than the dark matter halo radius of our galaxy M, then
we have that the companion galaxy is not close enough to significantly impact the interior mass,
and the interior mass is then simply given by the mass of the galaxy M.
2
def InteriorMass ( s e l f , r ) :
4 " " " I n t e r i o r Mass r e t u r n s t h e mass o f t h e o b j e c t ( s t a r / g a l a x y ) based on i t s
p o s i t i o n in r e l a t i o n to the c e n t e r o f the galaxy """
6
indices = r < s e l f . rthalo
8
10 i n t m a s s = np . f u l l ( r . shape , 0 . )
12 i f i n t m a s s [ i n d i c e s ] . shape != ( 0 , ) :
#I f r< s e l f . r t h a l o , then t h e i n t e r a c t i n g o b j e c t ( s t a r / g a l a x y ) i s c o n f i n e d
w i t h i n t h e dark matter h a l o r a d i u s o f t h e galaxy , and t h e i n t e r i o r mass
14 # i s c a l c u l a t e d as f o l l o w s
i n t m a s s [ i n d i c e s ] = ( s e l f . v h a l o ∗ ∗ 2 ) ∗ ( r [ i n d i c e s ] ∗ ∗ 3 ) / ( ( s e l f . a h a l o+r [
i n d i c e s ] ) ∗∗2)
16 i f i n t m a s s [ ~ i n d i c e s ] . shape != ( 0 , ) :
#Otherwise , t h e i n t e r i o r mass i s s i m p l y t h e mass o f t h e g a l a x y
18 intmass [~ i n d i c e s ] = s e l f . galmass
20 return intmass
Listing 7: Instance method "Interior Mass" used to calculate the mass of the interacting object
Next, we have an instance method called Density, which calculates the density of the surrounding
area of the galaxy (the density of the area it is moving through). This density is then used to
calculate the Dynamical Friction. The density instance method can be found in listing 8
def Density ( s e l f , r ) :
2 " " " Determines t h e d e n s i t y based on t h e i n p u t r which i s t h e p o s i t i o n o f
the o b j e c t in r e l a t i o n to the c e n t e r o f the galaxy . """
rinner = r ∗0.99
4 router = r ∗1.01
minner = s e l f . I n t e r i o r M a s s ( r ∗ 0 . 9 9 )
6 mouter = s e l f . I n t e r i o r M a s s ( r ∗ 1 . 0 1 )
dm = ( mouter−minner )
8 v o l =(4/3) ∗np . p i ∗ ( ( r o u t e r ∗ ∗ 3 ) −( r i n n e r ∗ ∗ 3 ) )
d e n s i t y=dm/ v o l
10
8
return density
Listing 8: Instance method "Density" used to calculate the density of the surrounding medium
being traversed
The last instance method we use in the Galaxy class is called DynFric (Listing 9). Here,
an estimate is calculated of dynamical friction, through use of a slightly altered version of the
Chandrasekhar formula given in equation 1. The friction force is then calculated using equation
9. In this equation, v represents the differential vector between the speeds of two objects, while
dv represents the magnitude of v. It is also important to note that in this altered version of the
Chandrasekhar formula, we use the following values G=1.0, and ln Λ = 3.0
4πG ln ΛM ρdv
f =− (9)
(1 + v)3
d e f DynFric ( s e l f , pmass , ppos , p v e l ) :
2 " " "Dynamic F r i c t i o n i s used when c o n s i d e r i n g f r i c t i o n i n our model . I t
c a l c u l a t e s t h e f r i c t i o n based on t h e mass p o s i t i o n and v e l o c i t y o f t h e o b j e c t .
The r e s u l t a n t f r i c t i o n c o n t r i b u t e s t o t h e o v e r a l l a c c e l e r a t i o n o f t h a t o b j e c t
"""
4 G=1.0
lnGamma=3.0
6 dv = p v e l − s e l f . g a l v e l
v = np . l i n a l g . norm ( dv )
8 dr = ppos − s e l f . g a l p o s
r = np . l i n a l g . norm ( dr )
10 galrho = s e l f . Density ( r )
f r i c m a g = 4 . 0 ∗ np . p i ∗G∗lnGamma∗ pmass ∗ g a l r h o ∗v /((1+ v ) ∗ ∗ 3 )
12 f r i c t i o n = (−dv/v ) ∗ f r i c m a g
14 return f r i c t i o n
Listing 9: Instance method "DynFric" used to estimate Dynamical Friction
In figure 2, we can see the effects of dynamical friction, and how it affects the results of the
simulation. In the figure, we see two instances of our simulation, one with dynamical friction
(b) and one without (a). (If viewed digitally, the still image is replaced with an animation).
These two images show the importance of dynamical friction in causing the galaxies to merge.
In (a), we see that the galaxies have extended tidal tails, but the two galaxies are not being
pulled together and will only be lightly disrupted. In (b), we see that the two galaxies are pulled
towards each other, and the galactic centers start moving towards each other, and the galaxies
merge into one larger galaxy.
9
(a) (b)
Figure 2: (a) Simulation at around 200 Myr without dynamical friction, (b) Simulation at around
200 Myr with dynamical friction enabled. In digital versions, click to play the videos.
3.2.2 StarGalaxy Class

Along with the Galaxy class, two other classes are also present in the code for our simulation.
The first o f w hich i s c alled S tarGalaxy. T he S tarGalaxy c lass i s a s ubclass o f o ur g alaxy class.
A subclass inherits all of the attributes and methods associated with the parent class, while it
can also have attributes and methods of its own. This can be seen by looking at listing 10. The
super() function is used to give access to the methods and attributes of the Galaxy class to our
new StarGalaxy class. We call all the attributes from the Galaxy class again, as can be seen in
our initialize statement, as well as several new ones, diskSize, galtheta, galphi, and n.
c l a s s StarGalaxy ( Galaxy ) :
2 " " "A sub−c l a s s o f galaxy , used t o combine i n f o r m a t i o n from t h e s t a r s and t h e
g a l a c t i c c e n t e r t o produce a g a l a x y with o r b i t i n g s t a r s . " " "
d e f __init__ ( s e l f , galmass , ahalo , vhalo , r t h a l o , g a l p o s , g a l v e l , d i s k S i z e , g a l t h e t a ,
galphi , n) :
4 s u p e r ( ) . __init__ ( galmass , ahalo , vhalo , r t h a l o , g a l p o s , g a l v e l )
s e l f . diskSize = diskSize
6 s e l f . galtheta = galtheta
s e l f . galphi = galphi
8 self .n = n
10 #D e f i n e t h e s t a r p o s i t i o n , v e l o c i t y , and a c c e l e r a t i o n i n such a way t h a t

th ey can be c a l l e d upon l a t e r more e a s i l y .
s e l f . s t a r p o s = np . f u l l ( ( 3 , s e l f . n ) , 0 . )
12 s e l f . s t a r v e l = np . f u l l ( ( 3 , s e l f . n ) , 0 . )
s e l f . s t a r a c c = np . f u l l ( ( 3 , s e l f . n ) , 0 . )
Listing 10: Creating the "StarGalaxy" class with inherited attributes
Most of these new attributes are parameters that can be altered by the user. Disksize is a fixed
value which can not be altered by the user, and represents the radius in which stars are present
within the galaxy. Galtheta and galphi give the orientation of the galaxy, more specifically, the
inclination and slew of the galaxy in the orbital plane. The effects of changing the inclination
and slew of a galaxy can be seen in figure 9. n represents the amount of stars being simulated.
The stars are then evenly split between the two galaxies such that each galaxy has 0.5n stars.
10
At the end of listing 10, we also create star templates, in which we create n stars, each with
the position velocity and acceleration set to 0. These will then be initialized in a later instance
method.
(a) (b) (c)
Figure 3: (a) Inclination and Slew of red galaxy are both 0◦ , (b) Inclination of red galaxy changed
to 90◦ , (c) Inclination and Slew of red galaxy are both 90◦
Now that we have gone over the new attributes defined in StarGalaxy, we can begin to look
at what new instance methods are defined in this class. The StarGalaxy class contains an
instance method called "MoveStars". The "MoveStars" method operates exactly the same as the
"MoveGalaxy" method, with the only difference being that we are replacing the galaxy’s position,
velocity and acceleration with the star’s position, velocity and acceleration. This instance method
can be found in listing 11.
d e f MoveStars ( s e l f , dtime ) :
2 " " " Function f o r moving a s t a r o v e r time by c a l c u l a t i n g t h e new p o s i t i o n
and v e l o c i t y " " "
n e w s t a r p o s = s e l f . s t a r p o s + s e l f . s t a r v e l ∗ dtime + 0 . 5 ∗ s e l f . s t a r a c c ∗ (
dtime ∗ ∗ 2 )
4 n e w s t a r v e l = s e l f . s t a r v e l + s e l f . s t a r a c c ∗ dtime
6 s e l f . starpos = newstarpos
s e l f . s t a rv e l = newstarvel
Listing 11: Defining a new "MoveStars" method for our StarGalaxy class
The next instance method present in our StarGalaxy class is the InitStars method. This method
initializes the position and velocity of the stars within the galaxy and can be found in listing 12.
def InitStars ( s e l f ) :
2 " " " I n i t S t a r s i n i t i a l i z e s t h e s t a r s i n t h e g a l a x y . I t g e n e r a t e s n number o f
s t a r s and g i v e s them random p o s i t i o n s w i t h i n t h e d i s k s i z e o f t h e g a l a x y . The
v e l o c i t i e s o f t h e s t a r s a r e a l s o c a l c u l a t e d based on t h e p o s i t i o n s a s w e l l a s
t h e i n p u t a n g l e s p h i and t h e t a " " "
c o s p h i = np . c o s ( s e l f . g a l p h i )
4 s i n p h i = np . s i n ( s e l f . g a l p h i )
c o s t h e t a = np . c o s ( s e l f . g a l t h e t a )
6 s i n t h e t a = np . s i n ( s e l f . g a l t h e t a )
f o r i in range ( s e l f . n) :
8 bad = True
11
w h i l e bad :
10 x t r y = s e l f . d i s k S i z e ∗ ( 1 . − 2 . ∗ np . random . random ( ) )
y t r y = s e l f . d i s k S i z e ∗ ( 1 . − 2 . ∗ np . random . random ( ) )
12 r t r y = np . s q r t ( x t r y ∗∗2+ y t r y ∗ ∗ 2 )
i f ( r t r y < s e l f . d i s k S i z e ) : bad = F a l s e
14
ztry = 0.0
16 xrot = xtry ∗ cosphi + ytry ∗ sinphi ∗ costheta + ztry ∗ sinphi ∗ sintheta
y r o t = −x t r y ∗ s i n p h i + y t r y ∗ c o s p h i ∗ c o s t h e t a + z t r y ∗ c o s p h i ∗ s i n t h e t a
18 z r o t = −y t r y ∗ s i n t h e t a + z t r y ∗ c o s t h e t a
r o t = np . a r r a y ( [ x r o t , y r o t , z r o t ] )
20 s e l f . s t a r p o s [ : , i ] = r o t + s e l f . g a l p o s . r e s h a p e ( −1)
22 v c i r c = np . s q r t ( s e l f . I n t e r i o r M a s s ( r t r y ) / r t r y )
24 v x t r y = −v c i r c ∗ y t r y / r t r y
vytry = v c i r c ∗ xtry / r t r y
26 vztry = 0.0
28 vxrot = vxtry ∗ cosphi + vytry ∗ s i n p h i ∗ costheta + vztry ∗ s i n p h i ∗ s i n t h e t a

v y r o t = −v x t r y ∗ s i n p h i + v y t r y ∗ c o s p h i ∗ c o s t h e t a + v z t r y ∗ c o s p h i ∗ s i n t h e t a
30 v z r o t = −v y t r y ∗ s i n t h e t a + v z t r y ∗ c o s t h e t a
32 v r o t = np . a r r a y ( [ vxrot , vyrot , v z r o t ] )
s e l f . s t a r v e l [ : , i ] = v r o t + s e l f . g a l v e l . r e s h a p e ( −1)
34 s e l f . s t a r a c c = np . f u l l ( ( 1 , 3 ) , 0 . )
Listing 12: Instance method "InitStars" used to initialize the positions and velocities of the stars
within each galaxy
An initial position of the star is first taken by taking the diskSize attribute and multiplying it
by a random number between -1 and 1. This is done for both the position in x as well as the
position in y. The position in z is always simply set to 0. From here, we calculate the magnitude
of our position vector. If the magnitude of the position vector of the star lies within the diskSize
attribute, then we take this as our initial position. Otherwise, if the magnitude of our position
vector lies outside our diskSize attribute, then the process is simply repeated until the star lies
within it. From here, we work with our galtheta and galphi attributes to further workout the
position of our stars depending on the orientation of our galaxy. The position of the star changes
dependent on the rotation undergone due to the inclination and slew angle. The expressions for
this can be found in equations 10 - 12.The derivation of these formulas can be found in Appendix
A.
xr = x cos φ + y sin φ cos θ + z sin φ sin θ (10)

yr = −x sin φ + y cos φ cos θ + z cos φ sin θ (11)
zr = −y sin θ + z cos θ (12)
After accounting for the rotations of our axes, we then simply add these positions to that of
the galaxy to obtain the positions of the stars within the galaxy. Next, we want to give these
stars an initial velocity. The first step to this is to calculate an estimate of the escape velocity
of the stars. This is done using an approximation of the escape velocity equation (equation 13).
In this equation, M is the interior mass which we calculate using the instance method defined
in the Galaxy Class. In our Galaxy class, we have also scaled all our equations such that we
12
take the gravitational constant to be equal to 1. Therefore, we can simply exclude it from our
equation. Furthermore, in the escape velocity equation, there is normally a factor 2 present in
the numerator, which we do not include in our estimation of the escape velocity.
r r
GM M
vesc = = (13)
r r
From here, we simply calculate velocities in x and y using equations 14 - 15. We assign the
velocity in the z direction to be equal to 0.
−vesc y
vx = p (14)
x2 + y 2
vesc x
vy = p (15)
x2 + y 2
From here, we use equations 10 - 12 to once again, account for the rotation of our axes, this
time replacing our position vector with the velocity vector. We then add our new velocity vector
of our star with the velocity vector of our galaxy to obtain an initial velocity vector for a star
within the galaxy. The initial acceleration of the stars is simply set to 0 for the time being.
The last new instance method present in our StarGalaxy class is another method called "Scale-
Mass" which can be found in listing 13. This new method simply multiplies the disk size of our
companion galaxy by the square root of the mass ratio between the companion galaxy and the
main galaxy. It also contains the super() function, and inherits all the attributes associated with
the ScaleMass function present in the Galaxy class.
d e f s c a l e M a s s ( s e l f , massFact ) :
2 " " "The f u n c t i o n s c a l e m a s s c a l c u l a t e s t h e d i s k s i z e o f a companion g a l a x y
based on t h e mass r a t i o o f t h e companion g a l a x y t o t h e main g a l a x y " " "
s e l f . d i s k S i z e = s e l f . d i s k S i z e ∗np . s q r t ( massFact )
4 s u p e r ( ) . s c a l e M a s s ( massFact )
Listing 13: Instance method "ScaleMass" used to calculate a new diskSize
3.2.3 Orbit Class

The final class that we define is the orbit class. The orbit class calculates the initial position
and velocity of the two galaxies based on the input masses of the galaxies. In our simulation,
a parabolic orbit is assumed. The attributes used in the Orbit class are energy, rp & tp, which
give the distances of closest approach, eccentricity, m1 & m2, which represent the masses of
the galaxies, and the bod1pos/bod2pos as well as bod1vel/bod2vel which represent the position
and velocities of each respective galaxy. For a parabolic orbit, we have that the energy of the
trajectory is 0, and the eccentricity is 1. We then calculate the initial positions and velocities of
the galaxies assuming a parabolic orbit as can be seen in listing 14.
c l a s s Orbit :
2 " " "The O r b i t c l a s s c a l c u l a t e s i n i t i a l p o s i t i o n and v e l o c i t y o f t h e two
g a l a x i e s based on a p a r a b o l i c o r b i t " " "
d e f __init__ ( s e l f , energy , rp , tp , e c c e n t r i c i t y , m1, m2, bod1pos , bod2pos , bod1vel ,
bod2vel ) :
13
Thijs Verkade 3.3 Creating the Graphical User Interface
4 self . energy = energy

self . rp = rp
6 self . tp = tp
self . eccentricity = eccentricity
8 self . m1 = m1
self . m2 = m2
10 self . bod1pos = bod1pos
self . bod2pos = bod2pos
12 self . bod1vel = bod1vel
self . bod2vel = bod2vel
14 self . initOrbit ()
16
18 def initOrbit ( s e l f ) :
" " " I n i t O r b i t i n i t i a l i z e s t h e o r b i t o f t h e two g a l a x i e s based on t h e i n p u t
masses o f t h e g a l a x i e s
20 This f u n c t i o n assumes a p a r a b o l i c o r b i t and r e t u r n s t h e p o s i t i o n and
v e l o c i t y o f both g a l a x i e s " " "
#P a r a b o l i c O r b i t
22 mu = s e l f . m1 + s e l f . m2
24 p = 2∗ s e l f . rp
nhat = np . s q r t (mu/ ( p ∗ ∗ 3 ) )
26 c o t s = 3 . 0 ∗ nhat ∗ s e l f . tp
s = np . a r c t a n ( 1 . 0 / c o t s )
28 c o t t h e t a = ( 1 . / ( np . tan ( s / 2 . ) ) ) ∗ ∗ 0 . 3 3 3 3
t h e t a = np . a r c t a n ( 1 . / c o t t h e t a )
30 t a n f o n 2 = 2 . / np . tan ( 2 . ∗ t h e t a )
r = ( p / 2 . ) ∗(1+ t a n f o n 2 ∗ ∗ 2 )
32
34 v e l = np . s q r t ( 2 . ∗mu/ r )
sinsqphi = p/(2.∗ r )
36 p h i = np . a r c s i n ( np . s q r t ( s i n s q p h i ) )
f = 2 . ∗ np . a r c t a n ( t a n f o n 2 )
38 xc = −r ∗np . c o s ( f )
yc = r ∗np . s i n ( f )
40 vxc = v e l ∗np . c o s ( f+p h i )
vyc = −v e l ∗np . s i n ( f+p h i )
42 xcom = s e l f . m2 ∗ xc / ( s e l f . m1+ s e l f . m2)
ycom = s e l f . m2 ∗ yc / ( s e l f . m1+ s e l f . m2)
44 vxcom = s e l f . m2 ∗ vxc / ( s e l f . m1 + s e l f . m2)
vycom = s e l f . m2 ∗ vyc / ( s e l f . m1+ s e l f . m2)
46
self . bod1pos = np . a r r a y ( [ [ − xcom ] , [ − ycom ] , [ 0 . 0 ] ] )
48 self . bod1vel = np . a r r a y ( [ [ − vxcom ] , [ − vycom ] , [ 0 . 0 ] ] )
self . bod2pos = np . a r r a y ( [ [ xc−xcom ] , [ yc−ycom ] , [ 0 . 0 ] ] )
50 self . bod2vel = np . a r r a y ( [ [ vxc−vxcom ] , [ vyc−vycom ] , [ 0 . 0 ] ] )
Listing 14: Initializing our "Orbit" class and the instance method "initOrbit" which calculates
the initial positions and velocities of two galaxies
3.3 Creating the Graphical User Interface

Now that we have established all our classes and methods that will simulate the galaxy interac-
tion, we must create a Graphical User Interface (GUI). The purpose of the GUI is to provide the
14
Thijs Verkade 3.3 Creating the Graphical User Interface
user an easy way to interact with the free parameters, as well as visually see how the galaxies
and stars interact with one another. Furthermore, the GUI is also used to provide important
information such as the time of the simulation, the separation of the galaxies at that time as well
as the relative velocity of the galaxies. Creating the User Interface is quite simple. With the help
of an application called Qt Designer, one can simply drag and drop widgets on a blank canvas,
and arrange them however they see fit. An blank canvas, along with the available widgets can
be found in figure 4.
Figure 4: A blank canvas in the Qt Designer Application
Using Qt Designer, to design the widgets, and then launching the application, we obtain the
GUI seen in figure 5
Figure 5: The designed GUI for the GalCrash Python application
15
Thijs Verkade 3.4 Combining the code with the GUI
3.4 Combining the code with the GUI

Having created the GUI, we now need to link the GUI with the code such that pushing a button,
or checking a box, actually has an impact. We need to write code which communicates between
our GUI and our classes such that we can send inputs from our GUI to use this information for
our parameters. The code should then calculate the relevant information about the stars and
galaxies, and then send this information back to the GUI so that it can be presented to the user
in an understandable way. To do this, we create a class called "GUI", which is created as a
subclass of "QMainWindow". The QMainWindow widget is a built-in class from PyQt which we
will use to link our widgets with our other classes.
c l a s s GUI( QMainWindow ) :
2 d e f __init__ ( s e l f ) :
s u p e r (GUI , s e l f ) . __init__ ( )
4 u i c . l o a d U i ( " GalCrashbeta . u i " , s e l f )
s e l f . show ( )
6 #When t h e i n i t i a l i z e button i s c l i c k e d , we c r e a t e t h e g a l a x y and s t a r s a s
w e l l a s g i v e an o r b i t u s i n g t h e f u n c t i o n s MakeGalaxy and MakeOrbit
s e l f . b u t t o n i n i t . c l i c k e d . c o n n e c t ( s e l f . MakeGalaxy )
8 s e l f . b u t t o n i n i t . c l i c k e d . c o n n e c t ( s e l f . MakeOrbit )
10
#When t h e s t a r t button i s p r e s s e d , we s t a r t s e e i n g an animated graph i n
which we s e e how t h e g a l a x y and s t a r s e v o l v e o v e r time and how they c o l l i d e
12 s e l f . b u t t o n s t a r t . c l i c k e d . c o n n e c t ( s e l f . update_plot )
14 #The r e s e t button c l e a r s t h e graph and new p a r a m e t e r s can be s e l e c t e d and

t h e s i m u l a t i o n can be run a g a i n
s e l f . restartbutton . c l i c k e d . connect ( s e l f . r e s e t )
16 s e l f . p a u s e b u t t o n . c l i c k e d . c o n n e c t ( s e l f . pause )
18 s e l f . time = 0 . 0 0
#C r e a t i n g l i s t s f o r which v a l u e s a r e appended t o when update p l o t i s
active
20 #These a r e then used t o c r e a t e t h e p l o t s o f t h e Galaxy S e p e r a t i o n and
R e l a t i v e V e l o c i t y o v e r Time
# C l e a r e d e v e r y time r e s e t i s p r e s s e d
22 s e l f . Dist = [ ]
s e l f . Vel = [ ]
24 s e l f . Time = [ ]
Listing 15: Creating a GUI class as a subclass of QMainWindow
3.4.1 Initializing the galaxies and stars

The first thing we do when creating our new GUI class is to load our GUI that we just built.
This will then load the general layout of our GUI as well as all the widgets we have created.
From here, we can start to send signals from the GUI to an instance method. Take for example,
our Initialize button, the code in lines 6-8 of listing 15 essentially say that when the initialize
button is pressed by the user, the "MakeGalaxy", as well as the "MakeOrbit" instance method
should be carried out. This is done for all the other buttons as can be seen in the following lines
of code.
Furthermore, we also set some initial variables at the top which are defined every time the
program is run. self.time is set to zero when the application is started and will be advanced in
16
the update plot function. self.Dist, self.Vel and self.Time are lists which will be filled with the
galaxy seperation, relative velocity of the galaxies and the time so that this information can be
plotted alongside the simulation.
Having connected all our buttons to instance methods, the next step is to set-up our graphs
as well as our text labels. This is just a matter of defining axis limits and labelling our axes.
The code for this can be found in listing 16. When this is done, our GUI looks ready to start
presenting information (see figure 6).
#g r a p h i c s V i e w i s t h e graph i n which t h e g a l a x i e s a r e d i s p l a y e d
2 s e l f . g r a p h i c s V i e w . setXRange ( −20 ,20)
s e l f . g r a p h i c s V i e w . setYRange ( −20 ,20)
4 s e l f . g r a p h i c s V i e w . h i d e A x i s ( ’ bottom ’ )
s e l f . graphicsView . hideAxis ( ’ l e f t ’ )
6
#v e l p l o t i s a p l o t o f t h e r e l a t i v e v e l o c i t i e s o v e r time
8 s e l f . v e l p l o t . setXRange ( 0 , 2 0 0 0 )
s e l f . v e l p l o t . setYRange ( 0 , 1 0 0 0 )
10 s e l f . v e l p l o t . s e t L a b e l ( ’ bottom ’ , ’ Time [ Myr ] ’ )
s e l f . v e l p l o t . s e t L a b e l ( ’ l e f t ’ , ’ R e l a t i v e V e l o c i t y [ km/ s ] ’ )
12
#C r e a t e d i f f e r e n t c o l o r s which a r e used f o r t h e a x e s
14 #mkPen o u t l i n e s t h e item i n t h e p r o p o s e d c o l o r
#mkBrush f i l l s t h e item i n t h e p r o p o s e d c o l o r
16
#Yellow
18 s e l f . pen1 = pg . mkPen ( c o l o r = ( 1 7 3 , 1 7 1 , 2 ) )
s e l f . brush1 = pg . mkBrush ( ’ y ’ )
20
#P u r p l e / Pink
22 s e l f . pen2 = pg . mkPen ( c o l o r = ( 1 5 5 , 4 1 , 1 6 3 ) )
24 #White
s e l f . pen3 = pg . mkPen ( c o l o r = ( 2 5 5 , 2 5 5 , 2 5 5 ) )
26
#Used f o r green text
28 s e l f . pen4 = pg . mkPen ( c o l o r = ( 8 5 , 1 7 0 , 0 ) )
#Used f o r red text
30 s e l f . pen5 = pg . mkPen ( c o l o r = ( 1 8 4 , 6 , 0 ) )
32 #l i g h t b l u e
s e l f . pen6 = pg . mkPen ( c o l o r = ( 2 8 , 1 7 6 , 2 1 7 ) )
34 s e l f . brush6 = pg . mkBrush ( c o l o r = ( 3 , 2 5 2 , 2 5 2 ) )
36 #Changing a x i s c o l o r s
s e l f . v e l p l o t . p l o t I t e m . g e t A x i s ( ’ l e f t ’ ) . s e t P e n ( s e l f . pen2 )
38 s e l f . v e l p l o t . p l o t I t e m . g e t A x i s ( ’ bottom ’ ) . s e t P e n ( s e l f . pen3 )
40
#d i s t p l o t i s a p l o t o f t h e g a l a x y s e p a r a t i o n o v e r time
42 s e l f . d i s t p l o t . setXRange ( 0 , 2 0 0 0 )
s e l f . d i s t p l o t . setYRange ( 0 , 2 0 0 )
44 #P l o t t e d above t h e v e l o c i t y p l o t s o we h i d e t h e x−a x i s a s t h i s i s common
between t h e two
s e l f . d i s t p l o t . h i d e A x i s ( ’ bottom ’ )
46 s e l f . d i s t p l o t . s e t L a b e l ( ’ l e f t ’ , ’ Galaxy S e p a r a t i o n [ kpc ] ’ )
s e l f . d i s t p l o t . p l o t I t e m . g e t A x i s ( ’ l e f t ’ ) . s e t P e n ( s e l f . pen4 )
48
17
#Add l a b e l s t o our g r a p h i c s V i e w p l o t which g i v e t h e Time , Galaxy

S e p a r a t i o n , and R e l a t i v e V e l o c i t y
50 s e l f . d i s t l a b e l . s e t S t y l e S h e e t ( " c o l o r : rgb ( 8 5 , 1 7 0 , 0 ) " )
s e l f . d i s t l a b e l . s e t T e x t ( " Galaxy S e p a r a t i o n : " )
52 s e l f . v e l l a b e l . s e t S t y l e S h e e t ( " c o l o r : rgb ( 1 5 5 , 4 1 , 1 6 3 ) " )
s e l f . v e l l a b e l . setText ( " Relative Velocity : " )
54 s e l f . t i m e l a b e l . s e t S t y l e S h e e t ( " c o l o r : white " )
s e l f . t i m e l a b e l . s e t T e x t ( "Time : " )
Listing 16: Setting up the graphs and defining text labels
Figure 6: An image of the GUI after defining axis parameters and labels
Now that we have set-up the GUI to look the way we want it, we need to examine the instance
methods that are called upon when certain buttons are pressed, to see exactly how these call
upon our main code to simulate the galaxy interactions. The first button that is pressed when
the user wishes to run the code after filling in all their desired parameters is the Initialize button.
Pressing this button calls upon two methods, the first of which is called "MakeGalaxy". The
purpose of MakeGalaxy, as the name suggests is two create two galaxies as objects part of our
StarGalaxy class. The code for this can be found in listing 17.
d e f MakeGalaxy ( s e l f ) :
2 " " "The f u n c t i o n MakeGalaxy c r e a t e s two g a l a x i e s based on t h e p a r a m e t e r s
p r o v i d e d i n t h e UI" " "
#g a l r e p r e s e n t s t h e main g a l a x y w h i l e comp r e p r e s e n t s t h e companion g a l a x y
4 #F i r s t some g e n e r a l p a r a m e t e r s taken from C h r i s Mihos ’ v e r s i o n
galmass = 4 . 8
6 ahalo = 0.1
vhalo = 1.0
8 rthalo = 5.0
g a l p o s = np . f u l l ( ( 3 , 1 ) , 0 )
10 g a l v e l = np . f u l l ( ( 3 , 1 ) , 0 )
diskSize = 2.5
12
#The f o l l o w i n g v a l u e s a r e g i v e n by t h e i n p u t s t h e u s e r has p r o v i d e d
14 g a l t h e t a = i n t ( s e l f . greenThetabox . v a l u e ( ) )
18
g a l p h i = i n t ( s e l f . greenPhibox . v a l u e ( ) )
16 comptheta = i n t ( s e l f . redThetabox . v a l u e ( ) )
compphi = i n t ( s e l f . r e d p h i b o x . v a l u e ( ) )
18 g a l n = i n t ( 0 . 5 ∗ i n t ( s e l f . starsnum . v a l u e ( ) ) )
compn = i n t ( 0 . 5 ∗ i n t ( s e l f . starsnum . v a l u e ( ) ) )
20 s e l f . g a l = StarGalaxy ( galmass , ahalo , vhalo , r t h a l o , g a l p o s , g a l v e l , d i s k S i z e ,
galtheta , galphi , galn )
s e l f . comp = StarGalaxy ( galmass , ahalo , vhalo , r t h a l o , g a l p o s , g a l v e l , d i s k S i z e ,
comptheta , compphi , compn )
22
#I f b i g h a l o s i s enabled , then our p a r a m e t e r s a r e s l i g h t l y d i f f e r e n t
24 i f s e l f . checkbighalo . isChecked ( ) :
s e l f . gal . rthalo = 20.0
26 s e l f . g a l . g a l m a s s = ( s e l f . g a l . v h a l o ∗∗2 ∗ s e l f . g a l . r t h a l o ∗ ∗ 3 ) / ( ( s e l f . g a l
. a h a l o+ s e l f . g a l . r t h a l o ) ∗ ∗ 2 )
else :
28 s e l f . g a l . galmass = 4 . 8
s e l f . gal . rthalo = 5.0
30
#Mass r a t i o i s p r o v i d e d by t h e u s e r
32 massrat = f l o a t ( s e l f . massratio . value ( ) )
34 #Using t h e mass r a t i o , we can s c a l e t h e p a r a m e t e r s o f t h e companion g a l a x y

according to the f u n c t i o n scaleMass
s e l f . comp . s c a l e M a s s ( m a s s r a t )
36
#I f b i g h a l o s i s enabled , then our p a r a m e t e r s a r e s l i g h t l y d i f f e r e n t
38 i f s e l f . checkbighalo . isChecked ( ) :
s e l f . comp . r t h a l o = s e l f . comp . r t h a l o ∗ 4 . 0
40 s e l f . comp . g a l m a s s = ( s e l f . comp . v h a l o ∗∗2 ∗ s e l f . comp . r t h a l o ∗ ∗ 3 ) / ( ( s e l f .
comp . a h a l o+ s e l f . comp . r t h a l o ) ∗ ∗ 2 )
Listing 17: Instance Method "MakeGalaxy" which creates the two galaxies used in the simulation
based on the parameters given by the user
As can be seen from listing 17, the first thing we do is to provide some default parameters.
Note that we also set the initial positions and velocities of the galaxies to 0. The actual positions
and velocities of the galaxies are calculated in our "MakeOrbit" method. To create the two
galaxies, we take the users input for the inclination and slew of both galaxies, as well as the
number of stars in each galaxy. From here, we create two galaxies with identical parameters, the
only difference being the angles provided.
Next, we check if the user has opted for big halos. If the big halos check box is ticked, then the
dark matter halos are four times bigger and more massive. The next step is then to scale the
mass of the companion galaxy, depending on what the user gave as a mass ratio. This is done
through our scaleMass method which we defined in our stargalaxy class. After we have scaled
the mass of our galaxy, we will then also apply the necessary changes to the dark matter halos
if the user has chosen for this. From here, we now have two galaxies with all our desired initial
parameters, except for their initial positions and velocities. This is where our next method,
"MakeOrbit", comes into play (see listing 18).
d e f MakeOrbit ( s e l f ) :
2 " " "Make o r b i t c r e a t e s t h e p a r a b o l i c o r b i t f o r both g a l a x i e s . I t then used
t h i s o r b i t t o d e t e r m i n e t h e i n i t i a l p o s i t i o n and v e l o c i t y o f t h e g a l a x y
c e n t e r s . F i n a l l y , i t i n i t i a l i z e t h e s t a r s i n s i d e both g a l a x i e s and then p l o t s
t h e g a l a x i e s a l o n g with t h e i r o r b i t i n g s t a r s " " "
19
4 energy = 0
eccentricity = 1
6 rperi = 3.0
8
#Parameter p r o v i d e d by u s e r ( p r o v i d e s t h e d i s t a n c e o f c l o s e s t approach )
10 t p e r i = f l o a t ( s e l f . peribox . value () )
12 #Returns t h e i n i t i a l p o s i t i o n and v e l o c i t i e s o f our g a l a x y c e n t e r s

s e l f . c r a s h O r b i t = O r b i t ( energy , r p e r i , t p e r i , e c c e n t r i c i t y , s e l f . g a l . galmass ,
s e l f . comp . galmass , s e l f . g a l . g a l p o s , s e l f . comp . g a l p o s , s e l f . g a l . g a l v e l , s e l f . comp .
galvel )
14
#S e t s t h e i n i t i a l p o s i t i o n and v e l o c i t y o f t h e g a l a x y c e n t e r s
16 s e l f . g a l . s e t P o s v e l ( s e l f . c r a s h O r b i t . bod1pos , s e l f . c r a s h O r b i t . b o d 1 v e l )
s e l f . comp . s e t P o s v e l ( s e l f . c r a s h O r b i t . bod2pos , s e l f . c r a s h O r b i t . b o d 2 v e l )
18
20 #I f b l u e c e n t e r e d i s checked , then we a l t e r t h e X and Y r a n g e o f t h e graph

s o t h a t t h e b l u e ( main ) g a l a x y
# i s always i n t h e c e n t e r o f our p l o t
22 i f s e l f . bluecenterbox . isChecked ( ) :
s e l f . g r a p h i c s V i e w . setXRange ( s e l f . g a l . g a l p o s [ 0 , 0 ] − 2 0 , s e l f . g a l . g a l p o s
[0 ,0]+20)
24 s e l f . g r a p h i c s V i e w . setYRange ( s e l f . g a l . g a l p o s [ 1 , 0 ] − 2 0 , s e l f . g a l . g a l p o s
[1 ,0]+20)
#The same p r i n c i p l e h o l d s f o r y e l l o w c e n t e r e d
26 i f s e l f . yellowcenterbox . isChecked ( ) :
s e l f . g r a p h i c s V i e w . setXRange ( s e l f . comp . g a l p o s [ 0 , 0 ] − 2 0 , s e l f . comp . g a l p o s
[0 ,0]+20)
28 s e l f . g r a p h i c s V i e w . setYRange ( s e l f . comp . g a l p o s [ 1 , 0 ] − 2 0 , s e l f . comp . g a l p o s
[1 ,0]+20)
#P l o t s t h e p o s i t i o n o f t h e g a l a c t i c c e n t e r f o r both g a l a x i e s
30 #We g i v e t h e g a l a x i e s s l i g h t l y d i f f e r e n t c o l o r s t o t h a t o f t h e s t a r s t o
d i s t i n g u i s h more e a s i l y between t h e two
Gal = pg . S c a t t e r P l o t I t e m ( brush = s e l f . brush6 , s i z e =6)
32 Gal . s e t D a t a ( [ s e l f . g a l . g a l p o s [ 0 , 0 ] ] , [ s e l f . g a l . g a l p o s [ 1 , 0 ] ] )
34 Comp = pg . S c a t t e r P l o t I t e m ( brush= s e l f . brush1 , s i z e =6)

Comp . s e t D a t a ( [ s e l f . comp . g a l p o s [ 0 , 0 ] ] , [ s e l f . comp . g a l p o s [ 1 , 0 ] ] )
36
#P l o t s t h e p o i n t s onto our p l o t
38 s e l f . g r a p h i c s V i e w . addItem ( Gal )
s e l f . g r a p h i c s V i e w . addItem (Comp)
40
42
44 #I n i t i a l i z e s t h e p a r a m e t e r s f o r t h e s t a r s i n both g a l a x i e s
s e l f . gal . InitStars ()
46 s e l f . comp . I n i t S t a r s ( )
48
50 #P l o t s t h e s u r r o u n d i n g s t a r s f o r both g a l a x i e s
s e l f . g r a p h i c s V i e w . p l o t ( ( s e l f . g a l . s t a r p o s [ 0 , : ] ) , ( s e l f . g a l . s t a r p o s [ 1 , : ] ) , pen
= None , symbol = " t 1 " , symbolPen = s e l f . pen6 , s y m b o l S i z e = 3 )
52 s e l f . g r a p h i c s V i e w . p l o t ( s e l f . comp . s t a r p o s [ 0 , : ] , s e l f . comp . s t a r p o s [ 1 , : ] , pen =
None , symbol = " t 1 " , symbolPen = s e l f . pen1 , s y m b o l S i z e = 3 , f i l l =True )
20
54 #d i s t i s t h e g a l a x y s e p a r a t i o n
d i s t = 3 . 5 ∗ np . l i n a l g . norm ( ( s e l f . g a l . g a l p o s − s e l f . comp . g a l p o s ) )
56 #We then p r i n t t h e g a l a x y s e p a r a t i o n on our graph
s e l f . d i s t l a b e l . s e t T e x t ( " Galaxy S e p a r a t i o n : { : . 2 f } kpc " . format ( d i s t ) )
58
#v e l i s t h e r e l a t i v e v e l o c i t y o f t h e two g a l a x i e s
60 v e l = 2 5 0 . ∗ np . l i n a l g . norm ( ( s e l f . g a l . g a l v e l − s e l f . comp . g a l v e l ) )
#We then p r i n t t h e r e l a t i v e v e l o c i t y on our graph
62 s e l f . v e l l a b e l . s e t T e x t ( " R e l a t i v e V e l o c i t y : { : . 2 f } km/ s " . format ( v e l ) )
64 #P r i n t i n i t a l time o f t=0
s e l f . t i m e l a b e l . s e t T e x t ( "Time : 0 Myr" )
Listing 18: Instance Method "MakeOrbit" which sets the initial positions and velocities of the
galaxies as well as creating stars around the galaxies
Similarly to the MakeGalaxy method, MakeOrbit starts with a few fixed default parameters,
along with the fact that for a parabolic orbit, the eccentricity is 1 while the energy is 0. From
here, the user provides the distance of closest approach and using this value, we create an object
in the class orbit which immediately initializes our orbit providing initial positions and velocities
for our galaxies which we then set to our galaxies using the setPosvel method.
From there, MakeOrbit deals with the plotting of our galaxy positions. The first thing that is
done is a check to see if the user has indicated that they wish the simulation to be centered on
either of the galaxies. If this is the case, then the axis limits are redefined such that the position
of the specified galaxy is always at the center of our simulation. Then, both galaxy positions
are plotted on our graph, and we initialize the stars for each galaxy using the InitStars method.
We then plot the positions of our stars for each galaxy in the corresponding color. In this way,
one can better understand how the evolved galaxy after the collision is comprised of a mixture
of stars from both initial galaxies. Then, the final thing we do in MakeOrbit is to calculate the
galaxy separation as well as the relative velocity of the galaxies, and then print these out in our
text labels.
Having now initialized our system, we obtain a graph which shows the two galaxies along with all
their stars, as well as giving us the separation between the galaxies and their relative velocities.
An image of our graph after our galaxies and stars have been initialized can be found in figure 7.
21
Figure 7: Graph after initialize button is pressed, showing the galaxies and their stars
3.4.2 Evolving the galaxies and stars over time

The next step, is to then evolve our system over time. When the Start button is pressed, our
method updateplot will be carried out, which will evolve our system over time, allowing us to
see how the galaxies interact. The updateplot method can be found in listing 19
d e f update_plot ( s e l f ) :
2 " " "The update p l o t f u n c t i o n r u n s t h e s i m u l a t i o n . I t c a l c u l a t e s t h e
a c c e l e r a t i o n o f t h e g a l a c t i c c e n t e r s f o r both g a l a x i e s a s w e l l a s t h e s t a r s
w i t h i n t h e g a l a x i e s . These a r e then used t o move t h e g a l a x y c e n t e r s a s w e l l a s
t h e s t a r s . The p l o t i s then updated t o show t h i s e v o l u t i o n . " " "
4 #We s e t t h e r e s t a r t b u t t o n and p a u s e b u t t o n t o True here , w h i l e t h e s e a r e

a c t i v e , t h e p l o t i s c o n s t a n t l y updated
#I f t h e pause button i s p r e s s e d , t h e update p l o t f u n c t i o n s t o p s r u n n i n g
but we can s t i l l s e e our p l o t s
6 #I f t h e r e s e t button i s p r e s s e d , then we c l e a r a l l our g r a p h s and r e s e t
t h e time back t o 0
s e l f . r e s t a r t b u t t o n . c l i c k e d = True
8 s e l f . p a u s e b u t t o n . c l i c k e d = True
10 #dtime taken from C h r i s Mihos ’ v e r s i o n , used t o e v o l v e t h e system o v e r

time
dtime = 0 . 0 4
12 i f s e l f . restartbutton . clicked :
while s e l f . pausebutton . c l i c k e d :
14 #C l e a r s t h e g r a p h s s o t h a t t h e new p o s i t i o n s can be p l o t t e d
without overlapping the previous p l o t
s e l f . graphicsView . c l e a r ( )
16 s e l f . distplot . clear ()
s e l f . velplot . clear ()
18
#C a l c u l a t i n g t h e a c c e l e r a t i o n u s i n g our a c c e l e r a t i o n f u n c t i o n
20 s e l f . g a l . g a l a c c = s e l f . comp . A c c e l e r a t i o n ( s e l f . g a l . g a l p o s )
s e l f . comp . g a l a c c = s e l f . g a l . A c c e l e r a t i o n ( s e l f . comp . g a l p o s )
22
22 #d i s t i s t h e r e l a t i v e d i s t a n c e between t h e two g a l a x i e s
d i s t = 3 . 5 ∗ np . l i n a l g . norm ( ( s e l f . g a l . g a l p o s − s e l f . comp . g a l p o s ) )
24
26 #I f Dynamic F r i c t i o n i s p r e s e n t ( g i v e n by u s e r ) , then a c c e l e r a t i o n
is calculated slightly differently
i f s e l f . c h e c k f r i c t i o n . isChecked ( ) :
28 s e l f . g a l . g a l a c c = s e l f . g a l . g a l a c c + s e l f . comp . DynFric ( s e l f . g a l
. InteriorMass ( d i s t /3.5) , s e l f . gal . galpos , s e l f . gal . g a l v e l )
s e l f . comp . g a l a c c = s e l f . comp . g a l a c c + s e l f . g a l . DynFric ( s e l f .
comp . I n t e r i o r M a s s ( d i s t / 3 . 5 ) , s e l f . comp . g a l p o s , s e l f . comp . g a l v e l )
30
comacc = ( ( s e l f . g a l . g a l m a s s ∗ s e l f . g a l . g a l a c c ) + ( s e l f . comp . g a l m a s s ∗
s e l f . comp . g a l a c c ) ) / ( s e l f . g a l . g a l m a s s + s e l f . comp . g a l m a s s )
32 s e l f . g a l . g a l a c c = s e l f . g a l . g a l a c c − comacc
s e l f . comp . g a l a c c = s e l f . comp . g a l a c c − comacc
34
s e l f . gal . staracc = s e l f . gal . Acceleration ( s e l f . gal . starpos ) + s e l f .
comp . A c c e l e r a t i o n ( s e l f . g a l . s t a r p o s )
36 s e l f . comp . s t a r a c c = s e l f . comp . A c c e l e r a t i o n ( s e l f . comp . s t a r p o s ) +
s e l f . g a l . A c c e l e r a t i o n ( s e l f . comp . s t a r p o s )
38 #A f t e r t h e a c c e l e r a t i o n i s c a l c u l a t e d , we then e v o l v e t h e system
u s i n g t h e f u n c t i o n s MoveGalaxy and MoveStars
40
s e l f . g a l . MoveGalaxy ( dtime )
42 s e l f . g a l . MoveStars ( dtime )
44 s e l f . comp . MoveGalaxy ( dtime )

s e l f . comp . MoveStars ( dtime )
46
#d i s t and v e l a r e t h e r e l a t i v e p o s i t i o n and v e l o c i t y between t h e
two g a l a x i e s
48 d i s t = 3 . 5 ∗ np . l i n a l g . norm ( ( s e l f . g a l . g a l p o s − s e l f . comp . g a l p o s ) )
v e l = 2 5 0 . ∗ np . l i n a l g . norm ( ( s e l f . g a l . g a l v e l − s e l f . comp . g a l v e l ) )
50
52 #Again , i f b l u e c e n t e r e d i s checked , then we want t h e graph t o

c o n s t a n t l y have t h e b l u e ( main ) g a l a x y a t t h e f o c u s , s o we c o n s t a n t l y update
t h e X and Y a x i s l i m i t s s o t h a t t h e g a l a x y a p p e a r s a t t h e c e n t e r
i f s e l f . bluecenterbox . isChecked ( ) :
54 s e l f . g r a p h i c s V i e w . setXRange ( s e l f . g a l . g a l p o s [ 0 , 0 ] − 2 0 , s e l f . g a l .
galpos [0 ,0]+20)
s e l f . g r a p h i c s V i e w . setYRange ( s e l f . g a l . g a l p o s [ 1 , 0 ] − 2 0 , s e l f . g a l .
galpos [1 ,0]+20)
56 #S i m i l a r l y f o r t h e y e l l o w g a l a x y
i f s e l f . yellowcenterbox . isChecked ( ) :
58 s e l f . g r a p h i c s V i e w . setXRange ( s e l f . comp . g a l p o s [ 0 , 0 ] − 2 0 , s e l f . comp
. galpos [0 ,0]+20)
s e l f . g r a p h i c s V i e w . setYRange ( s e l f . comp . g a l p o s [ 1 , 0 ] − 2 0 , s e l f . comp
. galpos [1 ,0]+20)
60
#P l o t t i n g t h e p o s i t i o n s o f t h e g a l a x i e s e x a c t l y t h e same a s b e f o r e
62 Gal = pg . S c a t t e r P l o t I t e m ( brush= s e l f . brush6 , symbol= ’ o ’ , s i z e =6)
Gal . s e t D a t a ( [ s e l f . g a l . g a l p o s [ 0 , 0 ] ] , [ s e l f . g a l . g a l p o s [ 1 , 0 ] ] )
64 Comp = pg . S c a t t e r P l o t I t e m ( brush= s e l f . brush1 , symbol= ’ o ’ , s i z e =6)
Comp . s e t D a t a ( [ s e l f . comp . g a l p o s [ 0 , 0 ] ] , [ s e l f . comp . g a l p o s [ 1 , 0 ] ] )
66 s e l f . g r a p h i c s V i e w . addItem ( Gal )
s e l f . g r a p h i c s V i e w . addItem (Comp)
23
68
#L i k e w i s e , we a l s o p l o t t h e p o s i t i o n s o f t h e s t a r s e x a c t l y t h e
same a s b e f o r e
70 s e l f . graphicsView . plot ( ( s e l f . gal . starpos [ 0 , : ] ) ,( s e l f . gal . starpos
[ 1 , : ] ) , pen = None , symbol = " t 1 " , symbolPen = s e l f . pen6 , s y m b o l S i z e = 3 )
s e l f . g r a p h i c s V i e w . p l o t ( ( s e l f . comp . s t a r p o s [ 0 , : ] ) , ( s e l f . comp . s t a r p o s
[ 1 , : ] ) , pen = None , symbol = " t 1 " , symbolPen = s e l f . pen1 , s y m b o l S i z e = 3 )
72
74 #Here we append our g a l a x y s e p a r a t i o n and r e l a t i v e v e l o c i t y t o t h e

arrays d efine d at the s t a r t
#These a r e then used t o p l o t t h e g a l a x y s e p a r a t i o n and r e l a t i v e
v e l o c i t y a s a f u n c t i o n o f time
76 s e l f . D i s t . append ( d i s t )
s e l f . Vel . append ( v e l )
78 #R e s u l t i s m u l t i p l i e d by 12 t o be i n Myr ( taken from C h r i s Mihos ’
version )
s e l f . Time . append ( s e l f . time ∗ 1 2 )
80
#Updates t h e l a b e l s with t h e new v a l u e s o f each parameter
82 s e l f . d i s t l a b e l . s e t T e x t ( " Galaxy S e p a r a t i o n : { : . 2 f } kpc " . format ( d i s t
))
s e l f . v e l l a b e l . s e t T e x t ( " R e l a t i v e V e l o c i t y : { : . 2 f } km/ s " . format ( v e l )
)
84 s e l f . t i m e l a b e l . s e t T e x t ( "Time : { : . 2 f } Myr" . format ( s e l f . time ∗ 1 2 ) )
86 #P l o t s t h e r e l a t i v e v e l o c i t y and g a l a x y s e p a r a t i o n
s e l f . v e l p l o t . p l o t ( s e l f . Time , s e l f . Vel , pen= s e l f . pen2 , width = 2 ,
name = " R e l a t i v e V e l o c i t y " )
88 s e l f . d i s t p l o t . p l o t ( s e l f . Time , s e l f . Dist , pen= s e l f . pen4 , width = 2 ,
name = " Galaxy S e p a r a t i o n " )
90
#p r o c e s s E v e n t s u p d a t e s t h e p l o t s
92 QApplication . processEvents ( )
94 #Advance t h e time
s e l f . time += dtime
96 else :
None
Listing 19: Instance Method "updateplot" which evolves the system in time
The first thing we do in updateplot is to set the status of our restart button and pause button
to True. While both of these booleans are true, our application will continuously evolve the
system over time. If either the pause button, or the reset button is pressed, then the updateplot
method will cease to run. The specifics of what the restart button and the pause button do, will
be discussed in the next section.
When the updateplot method is running, the first thing it does is it clears all our graphs. If
this is not done, then as the system evolves, the new positions of our galaxies and stars will be
plotted on an already exisiting plot of their old positions. By clearing the graphs each time,
we can essentially have the galaxy and stars "move" instead. The next thing our updateplot
function then does is to calculate the acceleration of both galaxies. This is done using our
Acceleration method from our Galaxy class which we defined earlier. If the user enables friction,
then the effect of dynamic friction will also be calculated using the method DynFric. The effect
24
of dynamic friction will then be added to our acceleration to find a new acceleration for our
galaxies. From here, we then find comacc, which represents the acceleration of the center of
mass of the two galaxies. This is then subtracted from each galaxies acceleration to obtain their
final acceleration.
Then, we calculate the acceleration of our stars by summing the acceleration due to both galaxies
separately. After having calculated the accelerations of both our galaxies as well as our stars, we
can then evolve the positons and velocities of these using the methods MoveGalaxy as well as
MoveStars. To do this, we use a fixed default parameter dtime to advance the system through
time. Now that our stars and galaxies have new positions and velocities, we calculate the galaxy
separation and relative velocity of the galaxies again.
We then plot the new positions of the galaxies and the stars, again, checking to see if the user
opted for the simulation to be centered around either galaxy and plotting accordingly. The galaxy
and star positions are plotted in exactly the same way as in our initialization stage. However,
here we also create two new plots, as we plot the galaxy separation over time, as well as the
relative velocity of the two galaxies over time. Then, we simply update the plots and advance
the time by dtime. This process is then repeated until the application is interrupted by the user.
3.4.3 Interrupting the simulation

Now we have a working simulation that plots the evolution of the galaxies and their stars over
time, giving the user information on the galaxy separation as well as the relative velocity of the
two galaxies. The simulation will run indefinitely unless the program is closed or interrupted
by the user. Two buttons have been added for the user to interrupt the program. The first
button is the reset button. When the reset button is pressed, as mentioned earlier, our boolean
will change from True to False, and our program will no longer continue to run the method
updateplot. Instead, another method is called upon, namely the reset method (see listing 20).
def reset ( s e l f ) :
2 " " " R e s e t c l e a r s a l l our v i s u a l p l o t s and i n f o r m a t i o n a n d r e s e t s t h e time t o
0 , s o t h a t new p a r a m e t e r s can be g i v e n and t h e s i m u l a t i o n can be run a g a i n " " "
4 self . restartbutton . clicked = False

self . pausebutton . c l i c k e d = False
6 self . graphicsView . c l e a r ( )
self . g r a p h i c s V i e w . setXRange ( −20 ,20)
8 self . g r a p h i c s V i e w . setYRange ( −20 ,20)
self . velplot . clear ()
10 self . distplot . clear ()
self . Dist = [ ]
12 self . Vel = [ ]
self . Time = [ ]
14 self . time = 0 . 0 0
self . d i s t l a b e l . s e t T e x t ( " Galaxy S e p e r a t i o n : " )
16 self . v e l l a b e l . setText ( " Relative Velocity : " )
self . t i m e l a b e l . s e t T e x t ( "Time : " )
Listing 20: Instance Method "reset" used to reset the simulation
Pressing the reset button thus interrupts our updateplot method, and it clears all our graphs,
resetting the time to 0. From here, the user can now change parameters and then initialize the
system again, and new galaxies with different properties can now be simulated. This allows the
25
user to easily study the effect of changing different parameters without having to restart the
application every time.
The other button present is the pause button. The pause button is similar to the reset button
in that pressing it interrupts the updateplot method. Pressing the pause button also calls upon
the pause method (listing 21)
d e f pause ( s e l f ) :
2 " " "We c r e a t e a s e p e r a t e f u n c t i o n f o r when pause i s c l i c k e d t h a t s i m p l y
s e t s t h e p a u s e b u t t o n . c l i c k e d b o o l e a n t o F a l s e This then means t h a t our p l o t
w i l l no l o n g e r be updated u n t i l s t a r t i s p r e s s e d once a g a i n " " "
s e l f . pausebutton . c l i c k e d = False
4 #When s t a r t i s p r e s s e d again , we re−run our update_plot f u n c t i o n i n which
we i m m e d i a t e l y s e t t h e p a u s e b u t t o n b o o l e a n t o True again , s o t h a t our
s i m u l a t i o n w i l l c o n t i n u e t o run from where i t l e f t o f f .
Listing 21: Instance Method "pause" used to stop the simulation temporarily
The difference between the reset and the pause button is that pressing the pause button does
not clear any of the previous information, so the graphs still remain intact allowing the user to
examine the evolution of the system at a specific time more easily. Pressing start again will then
simply resume the simulation and will continue to evolve the system.
3.4.4 Running the Application

Now that we have all the functionality desired in our program, the only thing left is to add some
code to start up the GUI when the python file is run. This is done in listing 22.
# h e r e we check i f t h e r e a l r e a d y e x i s t s an Q A p p l i c a t i o n . I f not we c r e a t e a new
one .
2 i f not Q A p p l i c a t i o n . i n s t a n c e ( ) :
app = Q A p p l i c a t i o n ( s y s . argv )
4 else :
app = Q A p p l i c a t i o n . i n s t a n c e ( )
6
window = GUI ( )
8
10
# here the a p p l i c a t i o n i s s t a r t e d
12
14 app . exec_ ( ) #f o r u s e an i n i n t e r a c t i v e s h e l l
Listing 22: Code to launch the GUI
26
Thijs Verkade 4. Evaluating the educational value of the Application
4 Evaluating the educational value of the Application

In order to evaluate the quality of the application as an educational tool, the galcrash Python tool
was tested by a tutorial group consisting of 5 students and 1 teaching assistant. The aim of the
tutorial was to "report on the simulations of galaxy interactions and mergers, and why more and
more powerful computers are required to conduct larger and more accurate simulations". In this
section, the positive and negative feedback received will be discussed, and will be used to come
up with further improvements to be made to the application in order to enrich its educational
quality. The full evaluation report can be found in Appendix B.
4.1 Discussing the feedback from tutorial session

In this session, students were able to observe different types of mergers by playing around with
the parameters. This allowed students to discover the characteristic results of major and minor
mergers, noting that at the end of a galaxy interaction, you are left with either a single merged
galaxy, two galaxies or two galaxies that are kind-of joined together. Students also discovered
that dynamical friction was vital to the galaxies merging, but were unaware of what exactly
dynamical friction was. Furthermore, they learned that the distribution of stars in galaxies are
impacted from a distance, even when the galaxies are not "touching" each other. Students also
observed that the simulation has limitations, and has very different parameters than real-life
galaxies, such as the simulated galaxies only having 250 stars, when in reality, the number of
stars present in a galaxy is on the order of billions. This is important to realize as the simulation
makes many assumptions and simplifications that are not necessarily correct, and therefore the
simulation should be used to help explain observations and gain a basic understanding of galaxy
interactions.
The galcrash application proved to be very fruitful in helping students understand the basics of
galaxy interactions. A large part of this was due to the fact that the code was written in Python.
As the students were somewhat familiar with Python, it allowed them to look into the code to
search for answers to questions that arose. Properties like how the initial velocities are calculated
can easily be found in the code, allowing the students to better understand the calculations and
mathematics going on within the code to produce the simulations. Looking at the code of such
an application also provided students with an example of how to run code with a user interface
in Python, which could be helpful for students looking to make similar simulations in the future.
Some criticisms of the application were made, one of which was that when the masses of the two
galaxies were chosen to be very different, things were "not obeying the laws of physics". This is
largely due to the problem of time-steps used in simulations. The properties of the stars and the
galaxies are constantly changing but in a simulation, we only change these properties at certain
time-steps, keeping them constant in between these time-steps. This leads to inaccuracies in
simulations especially at larger mass differences, in which the changes in the properties of the
stars and galaxies are much more volatile. Using a smaller time-step would help in improving the
accuracy of the simulation but would require many more calculations, slowing down the speed
of the application. While this could be done, this is not necessarily ideal as the purpose of the
application is to allow students to simulate galaxy interactions, and should be easily accessible,
not requiring supercomputers to be able to run such an application. The application is already
quite taxing, with students having had trouble running the simulations when the number of stars
was very high due to the speed or the overheating of computers.
27
Thijs Verkade 4.2 Changes made to the application
Having said this, this "loss of physics" could also be used as a learning opportunity to teach
students about using time-steps within simulations and thus give students more insight into the
drawbacks of creating simulations and where the limitations of such simulations lie.
4.2 Changes made to the application

Along with feedback mentioned in the evaluation report in Appendix B, feedback was also given
on certain technical issues, as well as requesting extra features that would enhance the experi-
ence of the students. Based on that feedback, several changes and additions were made to the
application to improve the overall quality of the application.
Previously, an initialize button was used to plot the galaxies based on the parameters provided
by the user. This was done is such a way that every time a user would change a parameter, they
would have to press Initialize again to see how changing that parameter reflects the initial state
of the simulation. This could make it quite difficult to see how certain parameters like inclination
and slew angles affected the orientation of the galaxies. Thus, the initialize button was removed,
and the application was reworked in such a way that the application automatically plots both
galaxies along with their stars based on the default parameters when launched. Then, when the
user provides any change in a parameter, this is also immediately reflected in the image of the
two galaxies, so the user has a more dynamic view of the impact of changing parameters on
the system. Furthermore, other small additions were also added to improve the quality of the
application. The application contains an option to have the graph centered around the main
galaxy. However, there was no option for the graph to be centered around the companion galaxy.
This is something students also wanted to do, and so this feature was added.
Another change made to the application was in the colour scheme. In the original design of
the application, the galaxies simulated were coloured red and green. However, upon receiving
feedback that such a scheme could prove quite difficult to follow for visually impaired or colour
blind people, the colour scheme was changed to a yellow-blue scheme instead. In this way, the
contrast was increased to improve overall visibility and make it more widely accessible to a variety
of students. There were also some minor technical issues with the GUI that were brought up
that are being reviewed.
4.3 Summary
Overall, the application proved to be a very useful tool for students to gain a basic understand-
ing of galaxy interactions and how changing certain parameters impacted the outcome of the
interaction. By having the application written in Python, it allowed students to be able to look
into the code to find more information about the underlying mathematics in such a manner that
it was understandable and could be followed.
Some concepts like major and minor mergers were described, as well as tidal tails. However,
students were unaware that they were describing these scenarios. Therefore, it may be useful
when using this application to introduce such concepts and vocabulary to the students, so that
they may better understand what is happening. Another important aspect that is not well
described is dynamical friction. As students discovered, dynamical friction is the driving force
behind the merging of galaxies. However, students were unaware of what dynamical friction
actually means. By having the concept of dynamical friction explained, students would be more
aware of why the galaxies are merging. Thus, the application would work best when accompanied
28
Thijs Verkade 4.3 Summary
with a tutorial, in which the concepts can be introduced and the students can explore these
concepts through help of the simulation.
29
Thijs Verkade 5. Conclusion
5 Conclusion
The aim of this project was to produce an updated version of the GalCrash application first
created by Chris Mihos to run using Java. The updated version was written in Python, with its
purpose being to educate students about galaxy interactions. The application runs simulations of
galaxy interactions between two spiral galaxies, taking parameter inputs from the user. The new
application written in Python has the same functionality as the old java application, but with
more explanations within the code, allowing the user to better understand how the underlying
code works, and what processes are taking place. The fact that the application is written in
Python also allowed students to look into the code to further explore the simulation while still
understanding what is being done.
According to the evaluation report, the application was overall very successful in introducing
students to the concepts of galaxy mergers. It allowed the students to explore parameters and
their effects on galaxy interactions. The application gave a very basic introduction to students
on the concepts of dynamical friction and allowed students to explore major and minor mergers.
However, the application lacked explanation as to what these concepts were. Students were able
to identify the impact of dynamical friction but were unaware of what it actually was. Students
discovered the properties of major and minor merger remnants without knowing what major and
minor mergers were. Thus while the application allowed students to explore galaxy interactions,
it did not offer a lot of information with regards to extra vocabulary and what the processes
taking place actually are. By accompanying the application with a tutorial, in which the concepts
are explained, the application would be much more successful in helping students understand
these concepts and then exploring the properties of these concepts.
5.1 Further Improvements

The new application does not add many new features that were missing in the java application.
This is something that could be improved upon further in future. Possible new features may
include things like presets, wherein one can easily choose two real galaxies to simulate their
merger. By selecting a preset, the user would not have to input any information about the
mass ratio or the number of stars, as this would already be known, but it would allow the user
to simulate the collision of two known galaxies easily. An example of this could be the Milky
Way-Andromeda galaxy collision. By having such presets, it allows the user to very easily select
their required collision parameters without having to do prior research on their properties.
Possible other features could see the inclusion of more diversity in the simulation. Being able to
choose different types of galaxies other than spirals, or different types of orbits for the galaxies
would allow the simulation much more flexibility in terms of its capabilities.
Another possible additional feature that could be added is to allow the user to specify a time
for which they want to see the evolution. As an example, if the user wishes to see how different
parameters affect the galaxy interaction at a time of 100Myr, then providing the user with such a
tool that they can specify the time for which they are interested in seeing the results. This could
make it easier for the user to compare the results of galaxy interactions for varying parameters.
30
Thijs Verkade A. Derivation of Rotation Matrix
A Derivation of Rotation Matrix

The GalCrash application allows for the user to input the inclination and slew of the two galaxies.
The inclination and slew of the galaxies affect the way in which the positions and velocities of
the stars within a galaxy are distributed. This is done via a rotation matrix. To better illustrate
our inclination and slew angles, we use figure 8.
Figure 8: Yaw, Pitch and Roll Rotation Axes (LaValle, 2006)
The yaw, pitch and roll axes all have well defined rotation matrices. By comparing our inclination
and slew axes to these three, we can determine the rotation matrices of our system. To do this,
we look at figure 9 again, which can also be found below.
(a) (b) (c)
Figure 9: (a) Inclination and Slew of red galaxy are both 0◦ , (b) Inclination of red galaxy changed
to 90◦ , (c) Inclination and Slew of red galaxy are both 90◦
Comparing these figures, we see that an increase in our inclination angle is a rotation around
the x-axis (roll) while an increase in our slew angle equates to a rotation around the z-axis (yaw).
However, yaw and roll occur when there is a counterclockwise rotation about the axes. In our
case, we are dealing with clockwise rotations, and so the sin terms in both our roll and yaw
matrix should be multiplied by -1. Our inclination angle is given by θ, and so, expressing the
31
Thijs Verkade B. Evaluation Report
roll axis rotation matrix in terms of θ, (and multiplying sin terms by -1), we obtain the following
matrix:
 
1 0 0
Ri = 0 cos θ sin θ  (16)
0 − sin θ cos θ
The slew angle is given in terms of φ. Writing the yaw rotation matrix in terms of φ (multi-
plying sin terms by -1), we have:
 
cos φ sin φ 0
Rs = − sin φ cos φ 0 (17)
0 0 1
   
cos φ sin φ 0 1 0 0
R = Rs × Ri = − sin φ cos φ 0 × 0 cos θ sin θ  (18)
0 0 1 0 − sin θ cos θ
Doing the multiplication, we obtain our rotation matrix:
 
cos φ cos θ sin φ sin θ sin φ
R = − sin φ cos θ cos φ sin θ cos φ (19)
0 − sin θ cos θ
To obtain the rotation in x,y, and z, we simply multiply by our position vector:
     
xr cos φ cos θ sin φ sin θ sin φ x
 yr  = − sin φ cos θ cos φ sin θ cos φ × y  (20)
zr 0 − sin θ cos θ z
Working this out then yields equations 10 - 12 used in the code (also seen below).
xr = x cos φ + y sin φ cos θ + z sin φ sin θ (21)

yr = −x sin φ + y cos φ cos θ + z cos φ sin θ (22)
zr = −y sin θ + z cos θ (23)
B Evaluation Report
32
Evaluation report – GalCrash Python
Situation:
The galcrash python tool was evaluated by a classroom observation and follow-up focus group.
The classroom situation was typical for a first-year introductory astronomy tutorial group at the
University of Groningen, though due to the COVID-19 restrictions in place at the time of the
study, the session was held online in Zoom. The 5 students and 1 Teaching Assistant (TA)
involved had all become familiar with using Zoom (or other online learning environments)
during the previous weeks of classes. The session was conducted in English, though no
participants were native English speakers, however, their regular language of instruction is
English, and all exhibit an English language level of at least C1 up to fluent. The students were all
at the end of their first year Bachelor studies in Astronomy at the University of Groningen. None
of the students, nor the TA, had participated in (or taught) previous versions of the tutorial
using pre-existing versions of the software. The session was independently observed by Dr Jake
Noel-Storr (reporting here), with the only contact during the session being to respond to
technical questions from the Teaching Assistant outside of the learning platform (as is the case
in the typical IA tutorial setup). The tutorial ran for 2 hours, and the focus group for 15 minutes
at the end.
The tutorial used is shown below:
33
Technical outcomes:
Technical issues and outcomes unrelated to learning were reported to the developer. None of
the extant technical issues had an observable effect on the overall learning experience – and
represented only recommendations for improvements to the user interface, or occurrences of
unexpected errors.
Observed outcomes:
The following learning outcomes were observed to take place during the classroom
observation. Paraphrased quotes from the session are included that highlight the learning that
took place, including from discussions between the students of things that they screenshotted
and shared as “interesting”.
Learning Outcomes; that… Exemplary Observations (paraphrased quotes)
Simmulations have limitations, but can help us to “When the masses are very different, things are
explain observations definitely not obeying the laws of physics”
“Wow – we are only using 250 stars, but galaxies

have billions”
34
Simmulations with more detail, require greater “Its so much slower when you go up to more
processing power and technology to complete stars”
them
“My laptop is melting now”
Astronomical images generally represent a single “Time is really slow”
snapshot in time
“We have to use the stop button to match an
image, so the image is not really the end point for
the galaxy interaction”
Astronomical images generally represent a single “Well we can move the simulation around to
two-dimensional representation of a three- make it look like the observation, but we cant
dimensional situation move the observation around so who knows in
3D if its right”
Looking at code can explain how a simulation “What are the initial velocities… oh wait we can
works and the assumptions that lie within it look at the code”
Galaxy mergers and interactions can have “We don’t know what dynamical friction does,
different outcomes but things won’t merge without it”
“At the end, sometimes they are oscillating (and

changing a bit each time), or sometimes a blob,
or sometimes are just carrying on forever”
Declared learning outcomes:

During the focus group, students were asked to describe what they thought that they had
learned, and these items are listed below:
⦁ That galaxies don’t directly crash into each other centre to centre
⦁ How to run code with a user interface in Python
⦁ Mass difference impacts the accuracy / reality (Note: this was made as a criticism of the
code, not the realization that is a problem with time-steps in simmulations)
⦁ At the end of an interaction or merger you can end up with just one thing, two things,
or two things that are kind-of joined together. Distances and velocities can become
zero, be oscillating, or keep moving away forever.
⦁ After mergers they can be joined together with spiral arms (Note: Describing Tidal
Tails – indicating that extra vocabulary should also be introduced first)
⦁ The distribution of stars in galaxies can be impacted even when the galaxies don’t
“touch” each other (i.e. Gravity acts at a distance).
⦁ That images are just 2-D, but each galaxy can be moved around in 3-D around 3 axes,
and then the view can also be moved in 3-D, so the physical situation is much more
complex than images reveal.
35
Thijs Verkade REFERENCES
References
Aceves, H., & Colosimo, M. (2006). Dynamical Friction in Stellar Systems: an introduction. arXiv
e-prints, arXiv physics/0603066, physics/0603066.
Delgado-Serrano, R., Hammer, F., Yang, Y. B., Puech, M., Flores, H., & Rodrigues, M. (2010).
How was the Hubble sequence 6 Gyr ago?, 509 arXiv 0906.2805, A78. https://doi.org/
10.1051/0004-6361/200912704
Dynamical friction. (n.d.). Swinburne University of Technology. https://astronomy.swin.edu.au/
cms/cpg15x/albums/userpics/dynamicalfriction1.gif
Fraknoi, A., Morrison, D., & Wolff, S. C. (2016). Galaxy mergers and active galactic nuclei
[Accessed = 12-06-2020]. https://openstax.org/books/astronomy/pages/28- 2- galaxy-
mergers-and-active-galactic-nuclei
Gallagher, J., & Sparke, L. S. (2007). Galaxies in the universe: An introduction (2nd ed.). Cam-
bridge University Press.
Lambas, D. G., Alonso, S., Mesa, V., & O’Mill, A. L. (2012). Galaxy interactions. I. Major
and minor mergers., 539 arXiv 1111.2291, A45. https://doi.org/10.1051/0004- 6361/
201117900
LaValle, S. M. (2006). Planning algorithms.
Mihos, C., Caley, D., & Bothun, G. (1999). Galcrash: N-body simulations on the student desktop.
36

bAST 2020 VerkadeT

Uploaded by

Copyright:

Available Formats

bAST 2020 VerkadeT

Uploaded by

Document Information

Original Title

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

bAST 2020 VerkadeT

Uploaded by

Copyright:

Available Formats

Rijksuniversiteit Groningen

Simulating Galaxy Collisions in Python

4 Evaluating the educational value of the Application 27

A Derivation of Rotation Matrix 31

(a) (b) (c)

3.2 Implementing OOP in Python

Listing 1: Creating a class with attributes

Listing 2: Creating Instance Methods to perform operations on our galaxy

Listing 5: Instance method "Acceleration" used to calculate acceleration of galaxy

2 def Potential ( s e l f , posin ) :

Listing 9: Instance method "DynFric" used to estimate Dynamical Friction

3.2.2 StarGalaxy Class

10 #D e f i n e t h e s t a r p o s i t i o n , v e l o c i t y , and a c c e l e r a t i o n i n such a way t h a t

Listing 10: Creating the "StarGalaxy" class with inherited attributes

(a) (b) (c)

28 vxrot = vxtry ∗ cosphi + vytry ∗ s i n p h i ∗ costheta + vztry ∗ s i n p h i ∗ s i n t h e t a

xr = x cos φ + y sin φ cos θ + z sin φ sin θ (10)

Listing 13: Instance method "ScaleMass" used to calculate a new diskSize

3.2.3 Orbit Class

4 self . energy = energy

3.3 Creating the Graphical User Interface

Figure 4: A blank canvas in the Qt Designer Application

Figure 5: The designed GUI for the GalCrash Python application

3.4 Combining the code with the GUI

14 #The r e s e t button c l e a r s t h e graph and new p a r a m e t e r s can be s e l e c t e d and

Listing 15: Creating a GUI class as a subclass of QMainWindow

3.4.1 Initializing the galaxies and stars

#Add l a b e l s t o our g r a p h i c s V i e w p l o t which g i v e t h e Time , Galaxy

Listing 16: Setting up the graphs and defining text labels

34 #Using t h e mass r a t i o , we can s c a l e t h e p a r a m e t e r s o f t h e companion g a l a x y

12 #Returns t h e i n i t i a l p o s i t i o n and v e l o c i t i e s o f our g a l a x y c e n t e r s

20 #I f b l u e c e n t e r e d i s checked , then we a l t e r t h e X and Y r a n g e o f t h e graph

34 Comp = pg . S c a t t e r P l o t I t e m ( brush= s e l f . brush1 , s i z e =6)

3.4.2 Evolving the galaxies and stars over time

4 #We s e t t h e r e s t a r t b u t t o n and p a u s e b u t t o n t o True here , w h i l e t h e s e a r e

10 #dtime taken from C h r i s Mihos ’ v e r s i o n , used t o e v o l v e t h e system o v e r

44 s e l f . comp . MoveGalaxy ( dtime )

52 #Again , i f b l u e c e n t e r e d i s checked , then we want t h e graph t o

74 #Here we append our g a l a x y s e p a r a t i o n and r e l a t i v e v e l o c i t y t o t h e

3.4.3 Interrupting the simulation

4 self . restartbutton . clicked = False

Listing 20: Instance Method "reset" used to reset the simulation

3.4.4 Running the Application

Listing 22: Code to launch the GUI

4 Evaluating the educational value of the Application

4.1 Discussing the feedback from tutorial session

4.2 Changes made to the application

5.1 Further Improvements

A Derivation of Rotation Matrix

Figure 8: Yaw, Pitch and Roll Rotation Axes (LaValle, 2006)

(a) (b) (c)

xr = x cos φ + y sin φ cos θ + z sin φ sin θ (21)

Evaluation report – GalCrash Python

“Wow – we are only using 250 stars, but galaxies

“At the end, sometimes they are oscillating (and

Declared learning outcomes:

⦁ How to run code with a user interface in Python

You might also like