Simple numerical model is developed to simulate magnetoelastic interactions in polycrystalline steel. Elastic properties and magnetic properties of ferromagnetic steel are related to each other, both of which depend on the magnetic domain structure. In the demagnetized state, six types of domains, each of which is magnetized parallel to one of the crystallographic axes, are randomly distributed in individual grain. The magnetic field increases the volume of the domains magnetized parallel to the field at the sacrifice of the other domains. When the domain realignment is almost completed, the domain magnetization starts to rotate. The stress also expands the domains parallel to the tension or perpendicular to the compression. In this study, we propose the simulation model which considers the volume change of the domains and the rotation magnetization, both induced by the magnetic field and the stress. The stress dependences of permeability and magnetostriction are predicted. Comparison of the results with the experiments ensures that the model is useful in explaining the experimental results.