Solution (source code)

= Solution

In the non-scattering <grey atmosphere> approximation, $I$ is the frequency-integrated <specific intensity> along a ray, $s$ is distance along that ray, $\rho$ is the <mass density>, and $\kappa$ is the <mass absorption coefficient>. The <local thermodynamic equilibrium> source is the frequency-integrated <Planck function> $B=\sigma_{\rm SB}T^4/\pi$. Across a thin element $ds$, absorption removes $\kappa\rho I\,ds$ and thermal emission adds $\kappa\rho B\,ds$. Their balance gives $dI/ds=-\kappa\rho(I-B)$.

Take the upward vertical as the polar axis, put $\mu=\cos\theta$, and define downward <optical depth> by
$$
\boxed{\tau(z)=\int_z^\infty\kappa(z')\rho(z')\,dz'.}
$$
Along a ray, $dz/ds=\mu$ and $d\tau/ds=-\kappa\rho\mu$. The <chain rule> therefore gives $\boxed{\mu I_\tau=I-B}$. In a plane-parallel, azimuthally symmetric field, the <solid angle> element is $d\omega=d\varphi\,d\mu$, so
$$
\boxed{J(\tau)=\frac12\int_{-1}^{1}I(\tau,\mu)\,d\mu.}
$$
<Local radiative equilibrium> requires zero net radiative heating of the material: absorption $4\pi\kappa\rho J$ balances emission $4\pi\kappa\rho B$. Where the absorption coefficient is nonzero this gives $B=J$, and substitution gives the required angular <radiative transfer equation>. It does not require the <specific intensity> itself to be isotropic.

At the top of an atmosphere with no external incoming long-wave radiation, impose $\boxed{I(0,\mu)=0\quad(-1\leq\mu<0)}$. These are inward directions; prescribing the same zero value for outward directions would remove the emergent radiation as well. The upward <radiative flux> and its derivative are
$$
\Phi=2\pi\int_{-1}^{1}\mu I\,d\mu,\qquad
\frac{d\Phi}{d\tau}=2\pi\int_{-1}^{1}(I-B)\,d\mu=4\pi(J-B)=0.
$$
Thus $\boxed{\Phi\text{ is constant with height}}$: <local radiative equilibrium> expresses energy conservation through successive atmospheric layers.