Added parameter "n_ray_pts" to distinguish it from "n_depth_pts" for clarity.

In the plane-parallel case the two are equal, but in, e.g., spherical
symmetry they are not.
Fixed full Lambda operator (almost).

Now when I operate Lstar (which is actually the full Lambda operator) on
the source function, it looks like regular Lambda iteration, except that
the depth boundary condition is wonky. Need to fix that...
Fixed some missing "elif"s
Changed "value" of the Planck function from 1 to 5.

This is to distinguish it from things which really do need to be one
(unity).
FINALLY FIXED APPROXIMATE LAMBDA OPERATOR!

The boundary value at depth was screwing up. It turns out for
upward-pointing rays, the coefficient for the matrix element at depth
needs to be 1, not 0. (This is confusing to me, given that Olson &
Kunasz (1987) says to "set the incident boundary intensities to zero".)

I'll need to think a lot more about why this is the right answer. But at
least for now it works. So I should be able to strip out any bandwidth
of this matrix to get a correct Lambda_star.