In this paper we investigate, through experiment and simulation, the effects of non-linear damping forces on the large amplitude structural dynamics of slender cantilever beams undergoing flapping motion in air. The aluminum beams are set into flapping motion through actuation at the beam base via a 4-bar crank-and-rocker mechanism. The beam strain response dynamics are investigated for two flapping amplitudes, 15° and 30°, and a range of flapping frequencies up to 1.3 times the first modal frequency. In addition to flapping at standard air pressure, flapping simulations and experiments are also performed at reduced air pressure (70% vacuum). In the simulations, linear and non-linear, internal and external damping force models in different functional forms are incorporated into a non-linear, inextensible beam theory. The external non-linear damping models are assumed to depend, parametrically, on ambient air density, beam width, and an empirically determined constant. Periodic solutions to the model equation are obtained numerically with a 1-mode Galerkin method and a high order time-spectral scheme. The effect of different damping forces on the stability of the computed periodic solutions is analyzed with the aid of Floquet theory. The strain-frequency response curves obtained with the various damping models suggest that, when compared to the linear viscous and non-linear internal damping models, the non-linear external damping models better represent the experimental damping forces in regions of primary and secondary resonances. In addition to providing improved correlation with experimental strain response amplitudes over the tested range of flapping frequencies, the non-linear (external) damping models yield stable periodic solutions for each flapping frequency which is consistent with the experimental observations. Changes in both the experimental ambient pressure and flapping amplitude are determined to result in some variation in the non-dimensional parameters associated with each of the non-linear external damping models. This result likely indicates an incomplete description of the model parameter dependence and/or non-linear functional form of the damping force.