# How do I interpret the energy per particle in the Smoothed Particle Hydrodynamics (SPH) formulation?

I'm working on a basic Smoothed Particle Hydrodynamics simulation. I figured I would start with an ideal gas with equation of state $$p=K\rho$$. However, in calculating the SPH equations of motion, we find $$\frac{\partial u}{\partial \rho}=p/\rho^2=K/\rho$$, so that $$u(\rho)=\text{const.}+K \log(\rho)$$, where $$u$$ is the energy per particle mass. This is a very weird energy equation, but if I plug it into my numerics I find that it does give conservation of energy! However I can't understand how to interpret this equation for $$u(\rho)$$.

Each particle ($$\bf a$$) in the simulation has mass $$m_a$$ and density $$\rho_a$$. We consider the volume as $$V_a=m_a/\rho_a$$, and so the thermodynamic change in internal energy is:

$$dU_a=T dS_a-p d\left(\frac{m_a}{\rho_a}\right)=T dS_a+\frac{p_a}{\rho_a^2} d\rho_a$$ We take $$m_a u_a=U_a$$ and $$dS=0$$ to get the formula for $$\partial u/\partial \rho$$. So I would naively expect that I could interpret the particle $$\bf a$$ to be a parcel of gas undergoing an isentropic process, but this gives a different dependence of $$u$$ on $$\rho$$!

It does make sense that my naive expectation doesn't have to be correct, because the particles in SPH don't have to be strict parcels of gas, they're merely points which are used to interpolate other fields like velocity and so on. But then, how do I interpret the energy per particle mass $$u_a(\rho_a)$$?

Extra information on the SPH Lagrangian:

Just to write down the specific system I'm talking about; If $$W_{ab}=W(r_a-r_b)$$ is the smoothing function, then I can write the Lagrangian for the system in its full form:

$$L=\sum_a \left(\frac{1}{2}m_a \dot{r}_a^2-m_a u_a(\sum_b m_b W_{ab})\right)$$

Where the second term is the energy per mass as a function of density, $$u_a(\rho_a)=C+K\log(\rho_a)$$ and the density is defined as $$\rho_ a=\sum_b m_b W_{ab}$$. This gives the SPH equations of motion:

$$\frac{DV_a}{Dt}=-\sum_b m_b\nabla_a W_{ab}\left(\frac{p_b}{\rho_b^2}+\frac{p_a}{\rho_a^2}\right)$$

with $$p_a=K\rho_a$$. So I can confirm in my simulations that with the SPH Lagrangian, this weird energy sum is conserved over time:

$$E=\sum_a \left(\frac{1}{2}m_a \dot{r}_a^2+m_a K\log(\sum_b m_b W_{ab})\right)$$

• By "this equation" you mean the energy one, $u(\rho)\sim\log(\rho)+\cdots$? Or one of the other two you've identified just beforehand? Commented Dec 5, 2023 at 12:50
• @KyleKanos thanks, fixed. My question was about the energy $u$. (in the context of SPH, the specific internal energy) Commented Dec 6, 2023 at 0:35

Answering my own question: That is correct for the equation of state $$p(\rho)=K\rho$$, and that is in fact the internal energy of a parcel of gas which I've forced to have that equation of state (EOS). But that is probably not the EOS that I want to use! Part of my confusion is that wasn't clear on what regime I'm interested in. All that we get from the SPH Lagrangian is that there's a function $$u$$ with units of energy per unit mass, and it's up to us to interpret the thermodynamics.
Say we consider an ideal gas at constant temperature. Then it's correct to use $$pV=nRT$$, so that $$p\propto \rho$$ and $$u(\rho)=C_1\log(\rho)+C_2$$. This formula is not as weird as I thought at first: it's the thermodynamic Helmholtz free energy $$A(T,V)=A(T,m/\rho)$$ per unit mass. The Helmholtz free energy is the appropriate thing to use when measuring mechanical work obtained from an isothermal system. Since we're now correctly interpreting this as isothermal, we can see that the entropy (as calculated from the Sackur-Tetrode equation) won't be conserved.
If we wanted to consider an ideal gas in the regime where thermal conduction is not an issue, we'd consider isentropic/adiabatic equation of state. The Sackur-Tetrode equation is a good starting point here, and section 43 of Landau and Lifshitz's Statistical Physics is a good reference. We get $$\frac{S}{k_B N}=\log(V)+c_V\log(U)+\text{stuff}$$. Setting $$dS=0$$ and using $$V=m/\rho$$ gives $$\frac{\partial U}{\partial \rho}=K \rho^{-1/3}$$. This is equal to $$p/\rho^2$$, and so we get equation of state $$p=K\rho^{5/3}=K\rho^\gamma$$. Here, $$u$$ corresponds to the thermodynamic internal energy $$U(S,V)=U(S,m/\rho)=C_1 \rho^{\gamma-1}+C_2$$. With this interpretation, the entropy is constant.