Monte Carlo (MC) simulations have been performed for the restricted primitive electrolyte model at a Bjerrum parameter B= 1.546 and three extremely small concentrations. The precision of the obtained thermodynamic properties (the excess energy, the heat capacity and the osmotic coefficient) is unprecedented in the literature. The objective is to study the first deviations from the Debye-Hückel limiting law (DHLL) without any assumptions. The values of the mentioned thermodynamic properties are significantly different from the values obtained from the extended Debye-Hückel expression (DH). MC values are situated between DH and DHLL values and are close to the values of the DHX theory (which uses a DH screened potential as the effective potential connected with the radial distribution functions). However, the MC values are also significantly different from the DHX values. The MC values of the radial distribution functions at contact are statistically identical to the values given by the DHX theory, but the contact values are quite uncertain because of the rare close encounters at extreme dilution. Nevertheless, the osmotic coefficients are well determined, since the contact term is vanishingly small in comparison to the third of the excess energy. The excess Helmholtz free energy has been calculated in the same Metropolis sampling by means of the Salsburg–Chesnut method. For an infinite Metropolis Markov chain, this method is shown to lead to the correct result for the excess electrostatic free energy. However, the method has not quite converged even with 10⁶ configurations, and the electrostatic free energies are found to be ca. 5% higher than the DHLL values. The DH theory predicts up to 5% lower values. Because of the present systematic deviation, the electrostatic entropies calculated by the difference between the excess energy and the electrostatic free energy are 15–20% lower than the DHLL values.