In order to understand the influence of magnetic fields on the propagation properties of waves, as derived from different local helioseismology techniques, forward modeling of waves is required. Such calculations need a model in magnetohydrostatic equilibrium as an initial atmosphere through which to propagate oscillations. We provide a method to construct such a model in equilibrium for a wide range of parameters, for use in simulations of artificial helioseismologic data. The method combines the advantages of self-similar solutions and current-distributed models. A set of models is developed by numerical integration of magnetohydrostatic equations from the subphotospheric to chromospheric layers.