In this paper large-eddy simulations (LES) of forced stratified turbulence using two common subgrid scale (SGS) models, the Kraichnan and Smagorinsky models, are studied. As found in previous studies using regular and hyper-viscosity, vorticity contours show elongated horizontal motions, which are layered in the vertical direction, along with intermittent Kelvin-Helmholtz (KH) instabilities. Increased stratification causes the layer thickness to collapse towards the dissipation scale, ultimately suppressing these instabilities. The vertical energy spectra are relatively flat out to a local maximum, which varies with the buoyancy frequency N. The horizontal energy spectra depend on the grid spacing Delta; if the resolution is fine enough, the horizontal spectrum shows an approximately -5/3 slope along with a bump at the buoyancy wavenumber k(b) = N/u(rms), where u(rms) is the root-mean-square (r.m.s.) velocity. Our results show that there is a critical value of the grid spacing Delta, below which dynamics of stratified turbulence are well-captured in LES. This critical Delta depends on the buoyancy scale L-b and varies with different SGS models: the Kraichnan model requires Delta < 0.47L(b), while the Smagorinsky model requires Delta < 0.17L(b). In other words, the Smagorinsky model is significantly more costly than the Kraichnan approach, as it requires three times the resolution to adequately capture stratified turbulence.