Monte Carlo (MC) simulations in ensembles with a fixed chemical potential or fugacity, for example the grand-canonical or the osmotic ensemble, are often used to compute phase equilibria. Chemical potentials can be computed either with an equation of state (EoS) or from molecular simulations. The accuracy of the computed chemical potentials depends on the quality of the (critical) parameters used in the EoS and the applied force field in the simulations. We investigated the consistency of both approaches for computing fugacities of the industrially relevant gases CO2, CH4, CO, H2, N2, and H2S. The critical temperature (Tc), pressure (Pc), and acentric factors (ω) of these gases are computed from MC simulations in the Gibbs ensemble. The effect of cutoff radius and tail corrections on the computed values of Tc, Pc, and ω is investigated. In addition, MC simulations in the Gibbs ensemble are used to compute the VLE of the 15 possible binary systems comprising the gases CO2, CH4, CO, H2, N2, and H2S, and the ternary systems CO2/CH4/H2S and CO2/CO/H2. Binary interaction parameters (kij) of these natural/synthesis gas mixtures are obtained by fitting the Peng-Robinson (PR) EoS to the binary VLE data from the MC simulations. The computed properties from the MC simulations are compared with the PR EoS, the GERG EoS, and experimental results. The MC results show that including tail corrections in the simulations is crucial to obtain accurate critical properties. The force fields used for the gases can reproduce the fugacities of the gases within 5% of the experimental data. The dew-point curves of all the 15 binaries were predicted correctly by the MC simulations, but the bubble-point curves for the systems H2/CO, CH4/H2, H2S/N2, and H2S/CO significantly deviate from the experiments.