  In neutrino oscillations, a neutrino created with one flavor can be later
detected with a different flavor, with some probability. In general, the
probability is computed exactly by diagonalizing the Hamiltonian operator that
describes the physical system and that drives the oscillations. Here we use an
alternative method developed by Ohlsson & Snellman to compute exact oscillation
probabilities, that bypasses diagonalization, and that produces expressions for
the probabilities that are straightforward to implement. The method employs
expansions of quantum operators in terms of SU(2) and SU(3) matrices. We
implement the method in the code NuOscProbExact, which we make publicly
available. It can be applied to any closed system of two or three neutrino
flavors described by an arbitrary time-independent Hamiltonian. This includes,
but is not limited to, oscillations in vacuum, in matter of constant density,
with non-standard matter interactions, and in a Lorentz-violating background.
