Solving differential equation for density in a gas sphere I'm trying to derive the equation for density as a function of height in a gas sphere due to gravitational force, and I have derived the following equation:
$$\frac{\mathrm{d}\rho}{\mathrm{d}z}=\frac{4 \pi G}{z^2RT}\rho(z)\int_{0}^{z}z'^2\rho(z')\mathrm{d}z'$$
Is there a way to solve it, even if numerically, for $\rho(z)$?
 A: My thought process for this is that I really wanted to get rid of that integral, so I started by isolating it on the right-hand side:
$$\dfrac{z^2RT}{4\pi G} \dfrac{\rho'}{\rho} = \int_0^z z'^2 \rho(z') dz'$$
Now we can differentiate both sides with respect to $z$, using the first part of the Fundamental Theorem of Calculus for the right-hand side.
$$\dfrac{2zRT + z^2RT'}{4\pi G} \dfrac{\rho'}{\rho} + \dfrac{z^2RT}{4\pi G} \dfrac{\rho\rho''-(\rho')^2}{\rho^2} = z^2\rho$$
$$(2zRT + z^2RT')\rho\rho' + z^2RT(\rho\rho''-(\rho')^2) = 4\pi Gz^2\rho^3$$
I definitely don't know how to solve that analytically, but with a given temperature profile you should be able to solve numerically using the form $\rho'' = \dfrac{4\pi G}{RT}\rho^2 + \dfrac{(\rho')^2}{\rho}-\dfrac{2T+zT'}{zT}\rho'$, from here I'd probably split it into a system and do Runge-Kutta, but the exact best way to do that is beyond my personal knowledge.
A: Here's a different rewrite of your integro-differential equation. I'm swallowing all of the constants up into one, called $C$.
$$\frac{1}{\rho(z)}\frac{d\rho}{dz} = \frac{C}{z^2}\int_0^z \zeta^2 \rho(\zeta)\,d\zeta $$
The left hand side can be seen as $\displaystyle\frac{d}{dz}\log(\rho(z))$. Differentiating both sides with respect to $z$ gives
$$ \frac{d^2}{dz^2}\log(\rho(z)) = -\frac{2C}{z^3}\int_0^z \zeta^2 \rho(\zeta)\,d\zeta + \frac{C}{z^2} z^2 \rho(z). $$
The middle term can be recognized as $\displaystyle-\frac{2}{z}\frac{d}{dz}\log(\rho)$, giving
$$ \frac{d^2}{dz^2} \log(\rho(z)) + \frac{2}{z} \frac{d}{dz}\log(\rho(z)) = C\rho(z). $$
This can be recast as
$$\frac{1}{z^2}\frac{d}{dz} \bigg(z^2\frac{d}{dz} \log(\rho(z))\bigg) = C\rho(z).$$
The operator on the left hand side is closely related to the radial Laplacian, which makes me think this might have been where you started? At any rate, this is horribly non-linear, so a numerical method is your best bet for sure.
A: Not a full answer, but here's what I have:
I isolated the integral, and then performed integration by parts twice.  I ended up defining a function $Q(z)$, such that $Q'''(z) = \rho(z)$.  This gives a nonlinear ODE:
$$\frac{z^2RT}{4\pi G}\frac{Q''''}{Q'''} = z^2Q'' - 2zQ' + 2Q(z)$$ where I used the arbitrariness of constants to define $Q(0) = 0$.  I can see that $\frac{d}{dz}\ln Q''' = \frac{Q''''}{Q'''}$, so if we went through the slog of integrating one more time, repeating the integration by parts trick, we would finally have an (I think) 5th order linear ODE with non constant coefficients.  I don't know how useful that is to you with the techniques you know, but my impression is that linear ODEs are far more tractable.
