The swimming of microorganisms typically involves the undulation or rotation of thin, filamentary objects in a fluid or other medium. Swimming in Newtonian fluids has been examined extensively, and only recently have investigations into microorganism swimming through non-Newtonian fluids and gels been explored. The equations that govern these more complex media are often nonlinear and require computational algorithms to study moderate to large amplitude motions of the swimmer. Here, we develop an immersed boundary method for handling fluid-structure interactions in a general two-phase medium, where one phase is a Newtonian fluid and the other phase is viscoelastic (e.g., a polymer melt or network). We use this algorithm to investigate the swimming of an undulating, filamentary swimmer in 2D (i.e., a sheet). A novel aspect of our method is that it allows one to specify how forces produced by the swimmer are distributed between the two phases of the fluid. The algorithm is validated by comparing theoretical predictions for small amplitude swimming in gels and viscoelastic fluids. We show how the swimming velocity depends on material parameters of the fluid and the interaction between the fluid and swimmer. In addition, we simulate the swimming ofCaenorhabditis elegans in viscoelastic fluids and find good agreement between the swimming speeds and fluid flows in our simulations and previous experimental measurements. These results suggest that our methodology provides an accurate means for exploring the physics of swimming through non-Newtonian fluids and gels.