Correlations measured in three dimensions in the Lyman-alpha forest are contaminated by the presence of the damping wings of high column density (HCD) absorbing systems of neutral hydrogen (H I; having column densities N(H I) > 1.6 x 10(17) atoms cm(-2)), which extend significantly beyond the redshift-space location of the absorber. We measure this effect as a function of the column density of the HCD absorbers and redshift by measuring threedimensional (3D) flux power spectra in cosmological hydrodynamical simulations from the Illustris project. Survey pipelines exclude regions containing the largest damping wings. We find that, even after this procedure, there is a scale-dependent correction to the 3D Lyman alpha forest flux power spectrum from residual contamination. We model this residual using a simple physical model of the HCD absorbers as linearly biased tracers of the matter density distribution, convolved with their Voigt profiles and integrated over the column density distribution function. We recommend the use of this model over existing models used in data analysis, which approximate the damping wings as top-hats and so miss shape information in the extended wings. The simple 'linear Voigt model' is statistically consistent with our simulation results for a mock residual contamination up to small scales (vertical bar k vertical bar < 1 h Mpc(-1)). It does not account for the effect of the highest column density absorbers on the smallest scales (e.g. vertical bar k vertical bar > 0.4 h Mpc(-1) for small damped Lyman-alpha absorbers; HCD absorbers with N(H-I) similar to 10(21) atoms cm(-2)). However, these systems are in any case preferentially removed from survey data. Our model is appropriate for an accurate analysis of the baryon acoustic oscillations feature. It is additionally essential for reconstructing the full shape of the 3D flux power spectrum.