The surfaces of neutron stars are likely sources of strongly polarized soft X rays due to the presence of strong magnetic fields. Scattering transport in the surface layers is critical to the determination of the emergent anisotropy of light intensity, and is strongly influenced by the complicated interplay between linear and circular polarization information. We have developed a magnetic Thomson scattering simulation to model the outer layers of fully-ionized atmospheres in such compact objects. Here we summarize emergent intensities and polarizations from extended atmospheric simulations, spanning considerable ranges of magnetic colatitudes. General relativistic propagation of light from the surface to infinity is fully included. The net polarization degrees are moderate and not very small when summing over a variety of field directions. These results provide an important foundation for observations of magnetars to be acquired by NASA's new IXPE X-ray polarimeter and future X-ray polarimetry missions.