matplotlib: changed to the more modern GridSpec method of rendering the plots
Removed line which prints J vector after every ALI iteration.
Reduced # of ALI iterations to 1 since we have full Lambda_star.

With full Lambda_star, doing ALI is equivalent to solving the integral version
of the radiative transfer equation, which requires no iteration but is
expensive. I did it this way just to make sure Lambda_star was being calculated
correctly.

Eventually we'll make the bandwidth of Lambda_star a knob. But as the code
stands currently it won't save any memory to reduce the bandwidth of
Lambda_star because we still call numpy.linalg.solve(), which assumes nothing
about the linear system while it's solving the system of equations. If we
wanted to save memory we would need to call, e.g., a tri-diagonal solver.
However I don't think Numpy comes with one built in.
Reduced the information in Lambda_star to just the tridiagonal elements.

This allows us to test whether our ALI scheme actually works. Since some
information is now missing from the Lambda operator, we need to iterate to find
global convergence.

I think we run into floating point issues before too long: after a few
iterations the innermost zones suddenly diverge to enormous values, even though
they behave perfectly well during the first few iterations. I noticed that
including more zones tends to delay the point at which the divergence occurs.
Replaced inner ray boundary condition with diffusion condition.

Currently this makes no difference in the solution since the atmosphere is
still isothermal. But it will make a big difference when the atmosphere is no
longer isothermal.