The system Eq. 27 for is non-linear with
respect to the and *n*_{e} because the coefficients of the
and -operators are quadratic in and the dependence of the Saha-Boltzmann factors and the collisional
rates on the electron density, respectively. The system is closed by the
abundance and charge conservation equations. To simplify the iteration
scheme, and to take advantage of the fact that not all levels strongly
influence all radiative transitions, we use a linearized and splitted
iteration scheme for the solution of Eq. 27. This scheme has
the further advantage that many different elements in different ionization
stages and even molecules can be treated consistently. A problem where
this is important is the modeling of nova and supernova atmospheres,
where there are typically very large temperature gradients within the
line forming region of the atmosphere.

To linearize Eq. 27, we follow [33] and replace
terms of the form
in Eq. 27 by
:

The most important advantage of this method is that it requires the solution of
large *linear* systems and low-dimensional non-linear system (for the electron
density). Thus, its solution is more stable and uses much less computer
resources (time and memory) than the direct solution of the original non-linear
equations. This allows us to treat many more levels with this method then
with more conventional algorithms. Using a nested iteration scheme like the
one described here will slow down the convergence of the iterations, but this
is more than offset for by the possibility of calculating much larger models
with less memory. Since we are able to solve a separate equation for each group
of elements, we can trivially parallelize the solution by distributing the
groups among the available processors.
Convergence acceleration methods can in principle be used, but they frequently
lead to convergence instabilities in the nested iterations for the solution of
the statistical equilibrium equations.

We have so far assumed that the electron density *n*_{e} is
given. Although this is a good assumption if only trace elements are
considered, the electron density may be sensitive to non-LTE
effects. This can be taken into account by using either a fixed point
iteration scheme for the electron density or, if many species or
molecules are included in the non-LTE equation of state, by a
modification of the LTE partition functions to include the effects of
non-LTE in the ionization equilibrium. The latter method replaces the
partition function, , with its non-LTE
generalization, , and uses
in the solution of the ionization/dissociation
equilibrium equation. We use this method because of the
large number of elements with various ionization stages as well
as molecules and condensation of dust grains included in
statistical equilibrium calculations (and not all of them in non-LTE).

Our iteration scheme for the solution of the multi-level non-LTE problem
can be summarized as follows: (1) for given *n*_{i} and *n*_{e}, solve the
radiative transfer equation at each wavelength point
and update the radiative rates and the
approximate rate operator, (2) solve the linear system Eq. 29
for each group for a given electron density, (3) compute new electron
densities (by either fixed point iteration or the generalized partition
function method), (4) if the electron density has not converged to
the prescribed accuracy, go back to step 2, otherwise go to step 1. The
iterations are repeated until a prescribed accuracy for the *n*_{e} and the
*n*_{i} is reached. It is important to account for coherent
scattering processes during the solution of the wavelength
dependent radiative transfer equation, it explicitly removes a global
coupling from the iterations.