pH is defined as the negative decadic logarithm of activity, a H = c H e βµ ex H , of hydronium ions, pH=-log 10 (a H /c ⊖ ), where c ⊖ = 1 M is the standard reference concentration, β = 1/k B T , and µ ex H is the excess electrochemical potential. The constant pH (cpH) Monte Carlo simulation method is a widely used approach for generating titration curves in systems undergoing protonation/deprotonation reactions . However, an indiscriminate application of this method poses a fundamental problem. In cpH simulations, entities such as proteins, colloidal particles, or polyelectrolytes are confined within a simulation box, while protons and ions have the freedom to exchange with an acid and salt in an implicit external reservoir . Consequently, the cpH simulation method is inherently semi-grand canonical. During the course of a cpH simulation, a proton is introa) Electronic mail: levin@if.ufrgs.br duced into the system from an external reservoir held at a predetermined pH. To maintain charge neutrality within the simulation cell, one of the cations or protons within the bulk of the cell is arbitrarily removed. However, this arbitrary removal lacks adherence to the principle of detailed balance, potentially yielding inaccurate outcomes , except in cases where the system contains a substantial amount of salt and is highly diluted in polyelectrolyte/protein. Fortunately, it is easy to rectify the standard cpH algorithm by incorporating a protonation step alongside a simultaneous grand canonical insertion of an anion. Conversely, a deprotonation step can be paired with a simultaneous grand canonical removal of an anion. This adjustment restores the detailed balance of the cpH algorithm, ensuring its internal consistency. We notice, however, that the simulation cell of such semi-grand canonical system will have a different mean electrostatic potential from that of the external reservoir. This is know as the Donnan potential. Therefore, the cpH simulation methods will allow us to predict the charge of polyelectrolyte only as a function of pH in the reservoir. Due to the presence of the Donnan potential, however, the pH of the reservoir, can be significantly different from that of an isolated system, for an identical protonation state . Thus, if one compares the titration curves, in which, say, the charge of colloidal particles is plotted as a function of the pH of the reservoir, with the titration curves of an isolated (canonical) system, there can be very large difference between the twoin particular for large volume fraction suspensions of low ionic strength . The difference between the two ensembles disappears in the limit or large ionic strength. This, perhaps, is the reason why this problem was not noticed previously -since most cpH simulations of proteins are performed at physiological concentrations and low protein volume fractions, when the difference between the two ensembles disappears. In fact, one can easily relate the pH c of a canonical system, in which the number of ions and protons is the same as the averages obtained using a semi-grand canonical simulation with a reservoir of pH gc , using equation