What do these negative mass values obtained from the Lane-Emdem equation entail? Starting from the hydrostatic equations for a star:
$$\frac{dM}{dr}=4\pi r^2\rho $$
$$\frac{dp}{dr}=-\frac{GM\rho}{r^2} $$
It is possible to obtain the following expression using the dimensionless variables from the Lane-Emdem equation:
$$\frac{M_n(\xi)}{M_n(\xi_1)}=\left(\frac{\xi}{\xi_1} \right)^2\frac{\phi_n(\xi)}{\phi_n(\xi_1)} $$
Where $ \xi_1$ is the "total radius" and $ \phi_n = d\theta_n/d\xi$. I've plotted this for $n=0,1,1.5,3,5$ and obtained:
The y axis says "Included mass/total mass" (left hand side of the equation). As you can see, some values are negative, for $n=1$ and $n=1.5$. This is because the $\phi_n$ is sometimes positive ($\phi_n(\xi_1)$ is negative in my case). But my question is, what does this represent? Are these, for example, prohibited regions where the star cannot have said radius? Or is it more likely to be some computational error? (By the way, I'm aware that the x axis should be labeled "xi". I mistook the letter earlier).
EDIT: Analytically, for $n=1$, $\theta(\xi)=\sin\xi/\xi$. Then:
$$\frac{M_n(\xi)}{M_n(\xi_1)}=\left(\frac{\xi}{\xi_1}\right)^2 \frac{1}{\phi_1(\xi_1)} \frac{\cos(\xi)\xi-\sin(\xi)}{\xi^2}$$
Which is not always positive. In other words, I don't think this is a numerical error.
 A: I asked one of my professors and the answer seems to be much more mundane than expected. Essentially, I was merely overshooting with my value for $\xi_1$. The value that one should use is simply such that $\theta(\xi_1)=0$; or in other words such that pressure is $0$ and no lower (negative). Without doing this, you are overextending the model into a non-physical region.
A: The Lane-Emden equation has, for some values of $n$, mathematical solutions where $\rho < 0$ in some regions of space. This then allows for the function $M(r)$ to have multiple zeroes and also negative values in some regions of space. This is possible because there is nothing in the equation constraining $\rho$ to non-negative values. So e.g. for $n=1$, we get an "onion" solution, with infinity of surfaces of zero density $\rho = 0$, where shells of positive density alternate with shells of negative density.
But in physics, we require that mass density must be positive or zero everywhere ($\rho \geq 0$). So whatever the solution that comes out of integration of the Lane-Emden equation, we believe it is relevant to models of stars only at those points where radius is between 0 and $r_1$, the first point where the density becomes zero. Then $r_1$ is called the "radius of the solution" and the rest of the solution for $r > r_1$ is pronounced unphysical, since we believe no material medium can have negative $\rho$.
Density outside the region $r\in \langle 0,r_1\rangle$ is not determined by the overall model, and it can be anything. Usually we think of an isolated sphere of gas in vacuum, in which case we put density outside to zero by hand.
