Use of amorphous silicon ($a$-Si) and hydrogenated amorphous silicon ($a$-Si:H) in photovoltaics has been limited by light-induced degradation (the Staebler-Wronski effect) and low hole mobilities, and voids have been implicated in both problems. Accurately modeling the void microstructure is critical to theoretically understanding the cause of these issues. Previous methods of modeling voids have involved removing atoms according to an {\it a priori} idea of void structure and/or using computationally expensive molecular dynamics. We propose a new fast and unbiased approach based on the established and efficient Wooten-Winer-Weaire (WWW) Monte Carlo method, by using a range of fixed densities to generate equilibrium structures of $a$-Si and $a$-Si:H that maintain 4-coordination. We find a smooth evolution in bond lengths, bond angles, and bond angle deviations $\Delta \theta$ as the density is changed around the equilibrium value of $4.9\times10^{22}\ $atoms/cm$^3$. However, a significant change occurs at densities below $4.3\times10^{22}\ $atoms/cm$^3$, where voids begin to form to relieve tensile stress, akin to a cavitation process in liquids. We find both small voids (radius $\sim$3 \AA) and larger ones (up to 7 \AA), which compare well with available experimental data. The voids have an influence on atomic structure up to 4 \AA beyond the void surface and are associated with decreasing structural order, measured by $\Delta\theta$. We also observe an increasing medium-range dihedral order with increasing density. Our method allows fast generation of statistical ensembles, resembles a physical process during experimental deposition, and provides a set of void structures for further studies of their effects on degradation, hole mobility, two-level systems, thermal transport, and elastic properties.