stellaGK
Making the parallel streaming in radially global stella implicit.
Integrating the exact quasineutrality solve with the response matrix approach is a work in progress that will soon be completely.
This however, doesn't include any radial variation from terms in the GKE proper. The naive way of including these (i.e. a big matrix) wouldn't work, since we really do rely on the speed of the tridiagonal solver for this term. What one could do is work in mixed x-ky space instead of pure Fourier. This won't work for any ky with connections, since this messes up the parallel boundary condition, but for periodic modes (e.g. zonal modes), this would be totally doable. The only catch is the gyroaveraging, which destroys the one-pointedness of the real-space representation. That can be worked around by using the un-gyroaveraged phi in the implicit routine (which then renders it one-point), and then adding the phi - < phi >_R term explicitly. This should remove the 1/k_perp^2 CFL constraint that's mentioned in the original Stella paper, at least for periodic modes.
