We calculate non-axisymmetric oscillations of neutron stars magnetized by purely poloidal magnetic fields. We use polytropes of index $n=1$ and 1.5 as a background model, where we ignore the equilibrium deformation due to the magnetic field. Since separation of variables is not possible for the oscillation of magnetized stars, we employ finite series expansions for the perturbations using spherical harmonic functions. Solving the oscillation equations as the boundary and eigenvalue problem, we find two kinds of discrete magnetic modes, that is, stable (oscillatory) magnetic modes and unstable (monotonically growing) magnetic modes. For isentropic models, the frequency or the growth rate of the magnetic modes is exactly proportional to $B_{\rm S}$, the strength of the field at the surface. The oscillation frequency and the growth rate are affected by the buoyant force in the interior, and the stable stratification tends to stabilize the unstable magnetic modes.

We calculate non-axisymmetric oscillations of uniformly rotating polytropes magnetized with a purely toroidal magnetic field, taking account of the effects of the deformation due to the magnetic field. As for rotation, we consider only the effects of Coriolis force on the oscillation modes, ignoring those of the centrifugal force, that is, of the rotational deformation of the star. Since separation of variables is not possible for the oscillation of rotating magnetized stars, we employ finite series expansions for the perturbations using spherical harmonic functions. We calculate magnetically modified normal modes such as $g$-, $f$-, $p$-, $r$-, and inertial modes. In the lowest order, the frequency shifts produced by the magnetic field scale with the square of the characteristic Alfv\'en frequency. As a measure of the effects of the magnetic field, we calculate the proportionality constant for the frequency shifts for various oscillation modes. We find that the effects of the deformation are significant for high frequency modes such as $f$- and $p$-modes but unimportant for low frequency modes such as $g$-, $r$-, and inertial modes.

We calculate axisymmetric toroidal modes of magnetized neutron stars with a solid crust in the general relativistic Cowling approximation. We assume that the interior of the star is threaded by a poloidal magnetic field, which is continuous at the surface with an outside dipole field. We examine the cases of the field strength $B_{{\rm{S}}}\sim10^{16}$ G at the surface. Since separation of variables is not possible for the oscillations of magnetized stars, we employ finite series expansions for the perturbations using spherical harmonic functions. We find discrete normal toroidal modes of odd parity, but no toroidal modes of even parity are found. The frequencies of the toroidal modes form distinct mode sequences and the frequency in a given mode sequence gradually decreases as the number of radial nodes of the eigenfunction increases. From the frequency spectra computed for neutron stars of different masses, we find that the frequency is almost exactly proportional to $B_{{\rm{S}}}$ and is well represented by a linear function of $R/M$ for a given $B_{{\rm{S}}}$, where $M$ and $R$ are the mass and radius of the star. The toroidal mode frequencies for $B_{{\rm{S}}}\sim 10^{15}$ G are in the frequency range of the quasi-periodic oscillations (QPOs) detected in the soft-gamma-ray repeaters, but we find that the toroidal normal modes cannot explain all the detected QPO frequencies.