  We calculate axisymmetric toroidal modes of magnetized neutron stars with asolid crust. We assume the interior of the star is threaded by a poloidalmagnetic field that is continuous at the surface with the outside dipole fieldwhose strength $B_p$ at the magnetic pole is $B_p\sim 10^{16}$G. Sinceseparation of variables is not possible for oscillations of magnetized stars,we employ finite series expansions of the perturbations using sphericalharmonic functions to represent the angular dependence of the oscillationmodes. For $B_p\sim 10^{16}$G, we find distinct mode sequences, in each ofwhich the oscillation frequency of the toroidal mode slowly increases as thenumber of radial nodes of the eigenfunction increases. The frequency spectrumof the toroidal modes for $B_p\sim 10^{16}$G is largely different from that ofthe crustal toroidal modes of the non-magnetized model, although the frequencyranges are overlapped each other. This suggests that an interpretation of theobserved QPOs based on the magnetic toroidal modes may be possible if the fieldstrength of the star is as strong as $B_p\sim 10^{16}$G.
We assume the interior of the star is threaded by a poloidal magnetic field that is continuous at the surface with the outside dipole field whose strength $B_p$ at the magnetic pole is $B_p\sim 10^{16}$G. Since separation of variables is not possible for oscillations of magnetized stars, we employ finite series expansions of the perturbations using spherical harmonic functions to represent the angular dependence of the oscillation modes.