The Peng-Robinson equation has been used to describe the phase behavior of elemental sulfur and the solubility of sulfur in natural gas mixtures. The sulfur is modelled as S-8 in all phases with b = 0.13122 m(3)/kmol and alpha(T) = 6.1051 + 2568.1/T MPa (m(3)/kmol)(2). Vapor pressures varying by 5 orders of magnitude are reproduced with a maximum error of 11% over a temperature range from 120 to 640 degrees C, and the calculated saturated liquid density is within 1% of the data up to 320 degrees C and within 3.6% at all temperatures. Calculated vapor densities are inaccurate because of the dissociation of Ss at low pressure or high temperature. The fugacity for pure solid sulfur is found from a separate model. A calculated phase diagram for sulfur-hydrogen sulfide matches three-phase lines and the solid-liquid-liquid-vapor quadruple point. Solubilities calculated from the proposed model show all the qualitative behaviors in the data and agree very well with experimental data at higher pressures in gases covering a broad range of compositions.