I have written my own gravitational engine! Although it started in python, I quickly ported it to C for a huge performance benefit.
The general idea is that it starts with a large array of 'stars'. Each is given some angular momentum. Then, the stars are gravitationally simulated over a period of time. This has led to some really cool effects!
Features:
Multithreading support
Customizable parameters through the command line
Autosaving and loading
Very low memory usage (Under 16MB for 100k stars)
High performance
Here is a simulation I did with n=1000 stars:
Here's another one with n=10k stars. This one turned out really well. The stars quickly form two galaxies, one much bigger than the other. Not only does the big one appear to be a spiral galaxy for a short period of time, but it also rips the smaller one apart!