I have been looking for a way to change some parameters included in the calculations of PySCF (e.g., changing the mass of the electron to the mass of a muon, so forming an exotic molecule). I know the whole code is on GitHub as it is a collaborative project, but I don't seem to find it. Has someone worked with something like this and knows where to find these parameters? Thanks!
-
1$\begingroup$ +1 and welcome to our new community! Thank you so much for contributing your question here and we hope to see much more of you in the future!!! Do you play chess by the way? $\endgroup$– Nike Dattani - No Free TimeCommented Nov 10, 2022 at 19:55
3 Answers
The electron mass is set in pyscf/data/nist.py along with other physical constants,
G_ELECTRON = 2.00231930436182 # http://physics.nist.gov/cgi-bin/cuu/Value?gem
E_MASS = 9.10938356e-31 # kg https://physics.nist.gov/cgi-bin/cuu/Value?me
Whether or not it's sufficient to change these is unclear to me. Note that the rest of the PySCF codebase appears to use the proton-to-electron mass ratio MP_ME = PROTON_MASS / E_MASS
rather than E_MASS
directly.
I'm afraid I don't think PySCF has the functionality... (I could be wrong!)
PySCF uses libcint to calculate all atomic orbital integrals. You would likely need to modify this codebase / find a different library that can calculate such integrals for muons. You could also write everything from scratch, but this would be quite involved.
Once you have the integrals you can feed them into PySCF and do all the chemistry stuff you'd like. As we normally work in atomic units (aka how the Hamiltonian is written), you will need to correct any terms where electron mass appears (replacing with that for muons). Terms not involving electrons, such as nuclear kinetic term will remain the same.
As already stated above, PySCF works with atomic units. This means that the electron mass is implicityly baked in the equations.
In the molecular Hamiltonian, the only term that depends on the mass is the kinetic energy. In atomic units, the kinetic energy is $\hat{T} = \hat{p}^2/2\mu = -\frac 1 2 \nabla^2$ where $\mu$ is the reduced mass $$ \mu = \frac {mM} {m+M}$$ which is chosen to be $\mu=1$ in atomic units. Note that for the hydrogen atom $m_e \ll m_p$ so the reduced mass is close to the electron mass $$ \mu = m_e \frac {1} {1 + \frac {m_e} {m_p}} = m_e \frac {1} {1 + 0.0005446170215} \approx 0.999455679425 m_e $$ but this won't be the case for heavier leptons as already for the muon one has $m_\mu/m_p \approx 0.1126095$.
To get the muonic operator, you would thus take $\hat{T} \to \mu_e/\mu_\mu \hat{T} $. Note, however, that the change of mass also significantly affects other parts of the calculation. The basis set requirements are different, since the muonic wave functions are more compact, and there will be significant muon-electron correlation effects.