Create Your Own Smoothed-Particle Hydrodynamics Simulation (With Python)

For today’s recreational coding exercise, we will simulate a star with smoothed-particle hydrodynamics (SPH). We will start with some initial condition and relax it into a stable stellar structure and measure the star’s density as a function of radius. The SPH method has many applications in astrophysics and elsewhere, and is in fact a general one to simulate all types of fluids.

You may find the accompanying Python code on github.

Before we begin, below is a gif of what running our simulation looks like:

SPH Method

We will represent the fluid as a collection of N point particles, indexed by i=1,…,N. Each particle has:

• mass m
• position rᵢ = [ xᵢ, yᵢ, zᵢ ]
• velocity vᵢ = [ vxᵢ, vyᵢ, vzᵢ ]

The particles will experience motion according to the Euler equation for an idea fluid:

Here, P is the fluid pressure and f is any additional external forces we may choose to impose on the system. We will assume that the pressure is given by a simple equation of state, called the polytropic equation of state (meaning it is just a function of density):