Chemical potential and relative chemical potential in grand canonical monte carlo

In the SemiGrandCanonicalEnsemble module, the documentation discusses about how the relative chemical potential is used to calculate potential difference and acceptance condition. However the input parameters of SemiGrandCanonicalEnsemble use chemical potential and judging by the code, this chemical potential is used the same way as relative chemical potential. So it is actually the relative chemical potential that is used during GCMC not chemical potential?

Hi @Jianli_Cheng,

The relative chemical potential is used. With two species A and B, you will get the exact same results if you use, for example
chemical_potentials = {'A': 0.5, 'B': 1.5}
as if you use
chemical_potentials = {'A': 10.5, 'B': 11.5}.

The acceptance criterion in the SGC ensemble is evaluated here.

@magnusrahm Thanks for the reply. However, I am still confused. I am running SGCMC to calculate the phase diagram of LiNiO2 as a function of temperature and Li concentration. I used a species X to represent Li vacancies and set chemical_potentials like this:
chemical_potentials = {'X': 0, 'Li': dmu}
Based on your answer, it seems the exact values of chemical_potentials don’t make. But later on, I will do free energy integration, which involves Li chemical potentials, and the free energy will change with different Li chemical potentials. So, how should I determine the input values for chemical_potentials?

I see, if the other species is vacancies and you put its chemical potential to zero, you could think of the provided chemical potential of Li as an absolute chemical potential. Note that the value of this chemical potentials is shifted by a constant if you have removed a slope from the energy as a function of Li concentration before fitting your cluster expansion (which one often does in order to have E=0 at c=0 and c=1, i.e., mixing energies).

I used the total energy calculated from DFT of each structure to fit my cluster expansion.
Based on your reply, I tested
chemical_potentials = {'X': 0, 'Li': -2}
chemical_potentials = {'X': -3, 'Li': -5}
X represents the Li vacancies. Indeed, I got the same results. Now, I want to calculate the grand canonical energy = E - mu_Li x. E is the potential (total energy) and x is the Li concentration. Both of them are calculated by SGCMC. So, what value should I use for mu_Li to calculate the grand canonical energy? Different mu_Li are clearly giving me different results. Should I calculate mu_Li by mu_Li=dE/dx?

Thank you so much for your time and reply!

I’m a little bit in deep waters here so take my words with a grain of salt. I think for most purposes it is just a choice of reference level that you could choose arbitrarily. However, if you want to relate your chemical potential to experiment conditions, you need to be more careful. Using \mu_\text{X}=0 is probably a good starting point, but you would also have to consider what your reservoir is. I did something similar for a system with hydrogen and vacancies, perhaps you could find some inspiration there (see especially Supplementary Note 4):

Thanks for the reply! I will think about it.

@magnusrahm After reading your paper, I am still not sure how to do free energy integration from SGCMC. Do I need to do VCSGC for free energy integration? It seems I cannot just use \mu_{Li} that was used as an input for SGCMC calculations.

I am following the method in this paper to do free energy integration. But I am not sure what \mu_{Li} I should use.

Hi @Jianli_Cheng

Integrating based on SGC simulations would be problematic if your system has two-phase region, such that one value of \mu corresponds to multiple concentrations, because your simulations will discontinously jump over the miscibility gap. VCSGC can solve that problem.

It is very hard to give an answer as to what \mu you should use because it depends on many details of how you have constructed your model and how you have done your simulations. I cannot give any specific advice there.