Describing neutrino oscillations in matter with Magnus expansion

We present new formalism for description of the neutrino oscillations in matter with varying density. The formalism is based on the Magnus expansion and has a virtue that the unitarity of the S-matrix is maintained in each order of perturbation theory. We show that the Magnus expansion provides better convergence of series: the restoration of unitarity leads to smaller deviations from the exact results especially in the regions of large transition probabilities. Various expansions are obtained depending on a basis of neutrino states and a way one splits the Hamiltonian into the self-commuting and non-commuting parts. In particular, we develop the Magnus expansion for the adiabatic perturbation theory which gives the best approximation. We apply the formalism to the neutrino oscillations in matter of the Earth and show that for the solar oscillation parameters the second order Magnus adiabatic expansion has better than 1% accuracy for all energies and trajectories. For the atmospheric $\Delta m^2$ and small 1-3 mixing the approximation works well ($< 3 %$ accuracy for $\sin^2 \theta_{13} = 0.01$) outside the resonance region (2.7 - 8) GeV.