Astrophysics and simulation
Simulation has always caught my eye. The exploration of systems outside of our own experience has always fascinated me. Oftentimes, though, it is difficult to see these other worlds with our own eyes, and numerical simulation on a computer can sometimes be the only approach to exploration. Simulation can be used for all sorts of purposes, but the simple visualization of complex phenomena is my favorite way to explore these other worlds. It can be visually stunning, as this project turned out to be.
I recently became interested in the simulation of astrophysics and astrophysical phenomena. I stumbled across a paper on the physics of gas flows, accretion disks, and other interesting physics. This began my study of fluid flows, and astrophysics.
This project site hosts an overview on my exploration of astrophysical phenomina and simulation. I’ve tried to make this visually interesting as well as intellectually, so feel free to just scroll down and look at the visuals.
During my various study, I came across PHANTOM, a free (!), open source (!), and highly accurate and configurable magnetohydrodynamic simulator. It uses smoothed particle hydrodynamics to accurately simulate a very wide range of astrophysical phenomena. This was a gold-mine for me, as it allowed me to do super accurate simulation with little barrier to entry. This is one of the more spectacular open-source finds I have come across.
After installing PHANTOM, I was ready to start simulating. PHANTOM is a smoothed particle hydrodynamics (SPH) code, widely used in published astrophysical research. The method is Lagrangian: rather than imposing a fixed grid and treating gas as a continuous field, PHANTOM discretizes the fluid into particles that carry the flow with them, each sampling local quantities such as density, pressure, and velocity. Because resolution naturally follows the mass, the code concentrates computational effort precisely where it is needed: collapsing cores, dense disk planes, shock fronts; rather than expending it on near-empty regions. This adaptivity, together with the conservation properties built into the SPH formulation, is what allows PHANTOM to resolve the full range of astrophysical behavior, from turbulence and shocks to the spiral structure of accretion disks and the evolution of magnetic fields.
As my first real custom simulation, I initialized a binary star system, surrounded by a large gaseous disk. This didn't take terribly long to simulate, which is surprising as it involved millions of particles and is simulated over the course of 618 years! It does this via some incredible maths and formidable calculations. It's incredible that one can simulate such astrophysical systems so accurately, and on a home computer nonetheless! For reference I performed this simulation, and most others on this page, on my Mac.
After it finished simulating, I needed a way to visualize the results. A simulation isn't much good if you can't analyze it. A common tool used in conjunction with PHANTOM is called SPLASH, a nifty tool to visualize smoothed-particle hydrodynamics. It's a little difficult to get going due to the limited documentation but once initiated it is very versatile. Most following image and video visualizations were generated using it.
Below is my simulation of the binary star system surrounded by a gaseous disk.
That was neat. After playing around a little more, I decided to try something larger. What I wanted to see was a tidal disruption event! This is when an object or more generally a fluid system engages in a very close flyby-event, where another massive object rips material away from the first object. This is called generally a tidal disruption, as tidal forces rip apart, or disrupt, the object. This is one of the more violent astrophysical phenomena.
Tidal disruption events, or TDE’s, can occur in a number of systems. The more violent of which occur in an actual star, where a star and its innards are ripped apart by another star or even black hole. As this TDE was one of my earlier experiments, and because star simulations are incredibly complicated, I decided to explore a protoplanetary disk disruption event.
A protoplanetary disk is just a bunch of matter and gas surrounding a star, like one might see in an early stage of a planetary system (proto - planetary). The disk, over a lengthy amount of time, would form planets and other, denser bodies surrounding the star, as the gas and matter coalesce. PHANTOM support a very rich variety of types of matter, where you can even simulate the dynamics of dust and asteroids among many, many other states. Again, PHANTOM is incredible.
I initialized two stars, both roughly the size of our sun on a course for a close fly-by event. This was a significant amount of compute for my poor computer. Fans were running full blast overnight, and it was sitting on icepacks I had in the freezer while running. In the mornings when I came to pause the simulation, every ice-pack was very thoroughly melted; they probably didn’t last more than a few hours.
This is a simulation of a fly-by event of 2 stars, one with a protoplanetary disk. The flyby should rip apart the protoplanetary disk, showing a clear example of tidal disruption from a close stellar flyby.
This simulation was run on my computer, on all cores. Every night for around a week straight I ran it when I had time, amounting to over 50 hours of simulation time, 48 frames were generated.
Simulated via PHANTOM, and visualized via splash. Files can be found below if you want to look at it youself.
In the beginning, you can see the star with the protoplanetary disk pulsate. This is merely stabilizing. if you look closely as it progresses, you may be able to see striations in the disk, forming natural spirals.
As we zoom in to the system just before the disruption event, you can see some fine structure present in the disk. In addition, notice the derformation of the disk due to the gravitation of the second star.
As the encounter progresses, the deformation becomes dramatic. The disk develops pronounced spiral arms and tidal tails. Material is clearly being ripped away, forming extended streams that trail off from the main disk. The central region remains bright, but the outer disk is now highly distorted and dynamically active. This formed very nice shapes. Also notice the turbulence around the second star, as matter and gas are violently flung around.
Now that the disruption has in large completed, around the original star you can see the remainder of the protoplanetary disk, with a significant amount of material ejection. In addition, around the second star, you can see a newly aquired and unstable disk of material ripped off from the primary star.
The same event at a different angle. This event occurs over the period of almost 7000 years!
Despite seeing these visualizations via splash in 3d, I still wanted to see them with a little more reality. It would be interesting to see if I could extract a 3d model from all my computer’s efforts. Maybe I could 3d print them as well.
Despite how they appear, these objects are quite nebulous. There is really no well defined edge or surface to these, due to their diffuse nature. The next best thing is to extract an isosurface from the cloud. An isosurface is a boundary enclosing or defining a surface of some constant value in a continuous field, and in this case an isosurface of constant density is probably what im looking for.
To do this I had to extract the data from the smoothed-particle representation of the mass distribution from PHANTOM. PHANTOM stores its simulated data in a series of “dump“ files and checkpoint files. These represent a similar format to video compression. The dump files represent keyframes, or a complete representation of all the data in the simulation. If the simulation stops, it can be restarted from these dump files. Similarly, the checkpoint files store the change in the data, how the particles moved and such; this is much less data-intensive.
For data extraction, the dump files are the easiest way to get the needed data. Fortunately, PHANTOM ships with a converter from smoothed-particle to a gridstream representation, which is a gridded array of values much more conducive to creating isosurfaces. After writing up some code (that you can find here, along with the rest of code for this project) to generate an isosurface, I chose the most visually interesting frame/dump file and created an isosurface enclosing some density around the primary star.
I also 3d printed it
Isosurface of primary star, at a lower density threshold
It is nice to see all that hard work done by my computer seen a little more concretely. Remember, this is a very accurate 7000 year long simulation, involving multiple millions of particles flying around stars and engaging in complicated fluid dynamics!
My next exploration was to explore the accretion disks around black holes. This would really push the limits of my computer, and even PHANTOM, in terms of accuracy. Here, effects are so strong that special relativity (and even general relativity!) play a significant role. Plasma surrounding black holes are subjected to unimaginable forces. Though PHANTOM does support relativistic effects, I wanted to use the most accurate simulation data to date for this.
To do this, I used data from simulations modeling the accretion disk around M87*—a supermassive black hole analogous to Sagittarius A* at the center of our own galaxy. The dataset, produced by the Illinois Relativity Group, is based on general relativistic magnetohydrodynamics (GRMHD), a framework that combines fluid dynamics, electromagnetism, and general relativity to model the behavior of magnetized plasma in very strong gravitational fields. This was perfect, as its also hosted online, with data free to download!
Here are some somewhat crude visualizations of the GRMHD v3 MAD-0.93 dataset, at varying times. Look at the detail!
Notice in red the areas of high density plasma falling in towards the black hole. Some interesting turbulence and magnetic effects there.
Also see in blue the low density regions at the poles of the black hole. This simulation managed to capture the relativistic jets coming off of the poles of the black hole! This is quite the feat!
Finally, I created a series of isosurfaces, 3d models, revealing some properties. Here, I have also included other interesting properties such as bulk velocity, magnetic feild strength, and more. There are many of these, but below are a few notable ones. Notice how the magnetic feild strength encloses a more globular region, yet with subtle trails. The isosurfaces with the spikes off the top and bottom are likely bulk momentum, as the plasma is getting shot out the relativistic jets at low density, yet imparted with enough speed to be comparable in momentum to the accretion matter just about the black hole!
iso-rho = 1.0E-01 • Black Hole
iso-rho = 1.5E-01 • Black Hole
iso-rho = 1.5E-01 • Black Hole
iso-u = 1.5E-05 • Black Hole
iso-u = 1.0E-06 • Black Hole
iso-rho_u = 1.0E-06 • Black Hole
I also 3d printed one of these. Turned out great!
I hope you enjoyed this exploration into astrophysics, I know I did. I think my poor computer needs a break from this sort of calculation for a while after what its been through. Overall, the simulations came out great.
