QMCPACK
Cutoff value in H.ccECP.xml
I think there is a problem in the file: https://pseudopotentiallibrary.org/recipes/H/ccECP/H.ccECP.xml In particular, the cutoff is set to: cutoff="1e-10". In fact, if I obtain the xml file from the UPF using ppconvert I obtain a cutoff of cutoff="0.807074".
By the way, it is not clear to me what is the meaning of this cutoff: is it the cutoff of the channel that is specified (I see there is a cutoff for s and a cutoff for p)? If so, the meaning of cutoff="1e-10" is that qmcpack do not see the pseudo of H and it is like having the Coulomb potential in H, is that correct?
@aannabe do you want to comment on the strategy here? All the other files created around the same time have the "expected" cutoff of O(1).
The cutoff values given in the *.xml files are not used in qmcpack so I don't think it matters and it shouldn't become just Coulomb. For all ccECPs, the cutoffs are radii where Coulomb starts to deviate within that radius from the semi-local potential for some threshold. But I am not sure yet why it's so small for H and He. I'll take a closer look and we will replace these if we find any pathologies.
@aannabe I am a little puzzled to know qmcpack ignores the cutoff. Why are them written in the xml file then? Also, this is from the qmcpack output:
Adding pseudopotential for H
Linear grid ri=0 rf=10 npts = 10001
ECPComponentBuilder::buildSemiLocalAndLocal
WARNING Nrule was not determined from qmcpack input or pseudopotential file. Setting sensible default.
Assuming Hartree unit
Number of angular momentum channels 2
Maximum angular momentum channel 1
Creating a Linear Grid Rmax=1e-10
Using global grid with delta = 0.001
Making L=1 a local potential with a radial cutoff of 9.999
Quadrature Nrule: 4
Non-local pseudopotential parameters
Maximum angular mementum = 0
Number of non-local channels = 1
l(0)=0
Cutoff radius = 0.001
Spherical grids and weights:
1 0 0 0.08333333333
-1 1.224646799e-16 0 0.08333333333
0.4472135955 0.894427191 0 0.08333333333
-0.4472135955 0.7236067977 0.5257311121 0.08333333333
0.4472135955 0.2763932023 0.8506508084 0.08333333333
-0.4472135955 -0.2763932023 0.8506508084 0.08333333333
0.4472135955 -0.7236067977 0.5257311121 0.08333333333
-0.4472135955 -0.894427191 1.095357397e-16 0.08333333333
0.4472135955 -0.7236067977 -0.5257311121 0.08333333333
-0.4472135955 -0.2763932023 -0.8506508084 0.08333333333
0.4472135955 0.2763932023 -0.8506508084 0.08333333333
-0.4472135955 0.7236067977 -0.5257311121 0.08333333333
Maximum cutoff radius 0.001
QMCHamiltonian::addOperator LocalECP to H, physical Hamiltonian
See that there is a line telling "Creating a Linear Grid Rmax=1e-10". What is the meaning of this line? Also notice the "Maximum cutoff radius 0.001".
In contrast, when I use the xml obtained from upf, where the cutoff is set to 0.807074, I have in output:
Adding pseudopotential for H
Linear grid ri=0 rf=10 npts = 10001
ECPComponentBuilder::buildSemiLocalAndLocal
WARNING Nrule was not determined from qmcpack input or pseudopotential file. Setting sensible default.
Assuming Hartree unit
Number of angular momentum channels 2
Maximum angular momentum channel 1
Creating a Linear Grid Rmax=0.807074
Using global grid with delta = 0.001
Making L=0 a local potential with a radial cutoff of 9.999
Quadrature Nrule: 4
Non-local pseudopotential parameters
Maximum angular mementum = 1
Number of non-local channels = 1
l(0)=1
Cutoff radius = 0.808
Spherical grids and weights:
1 0 0 0.08333333333
-1 1.224646799e-16 0 0.08333333333
0.4472135955 0.894427191 0 0.08333333333
-0.4472135955 0.7236067977 0.5257311121 0.08333333333
0.4472135955 0.2763932023 0.8506508084 0.08333333333
-0.4472135955 -0.2763932023 0.8506508084 0.08333333333
0.4472135955 -0.7236067977 0.5257311121 0.08333333333
-0.4472135955 -0.894427191 1.095357397e-16 0.08333333333
0.4472135955 -0.7236067977 -0.5257311121 0.08333333333
-0.4472135955 -0.2763932023 -0.8506508084 0.08333333333
0.4472135955 0.2763932023 -0.8506508084 0.08333333333
-0.4472135955 0.7236067977 -0.5257311121 0.08333333333
Maximum cutoff radius 0.808
Hi @zenandrea, there does seem to be an issue, thanks for pointing it out. I was looking at the source code and the potential files. Inside qmcpack, we actually do actually use the cutoff by finding the maximum cutoff across all potential channels, and then use that to define the grid for which the splined potentials are defined on for the non-local channels. This is supposed to be a purely local potential, where we only smooth out the coulomb potential at the origin to remove the divergence. However, something was overlooked when creating these xml files that is ultimately the source of your problems.
The reason things are off is the following. Both of these xml files were generated using ppconvert, which reads in the gaussian representation of the semilocal potential in a GAMESS format. If you notice, He for example has the format of
He-ccECP GEN 0 1 3 2.000000 1 32.000000 64.00000 3 32.000000 -27.70084 2 33.713355 1 0.000000 2 1.0000000
which has a 3 gaussian local channel to smooth out the coulomb potential and we add a zeroed out nonlocal S channel in order to get it to actually run in GAMESS/ppconvert/other codes (some don't allow you to just have a local channel for some reason, so we typically do this). However, when ppconvert encounters this, it will also include the zeroed out S channel which we don't want. We have to remove it after the fact to make sure we only have 1 local channel only in the xml file and no nonlocal channels.
ppconvert wrote the xml with potentials for S and P
<semilocal units="hartree" format="r*V" npots-down="2" npots-up="0" l-local="1"> and sets the local potential as P.
What we actually want for these is a purely local potential
<semilocal units="hartree" format="r*V" npots-down="1" npots-up="0" l-local="0"> and there should only be one vps xml tag instead of two. This was overlooked. Notice that your qmcpack output is including a quadrature rule for a potential that is supposed to purely local, which means the code is also allocating the nonlocal channel when we don't want it.
Basically, the xml generated potentials in order to get ppconvert to actually give you the right thing had to modified by hand after the fact to make sure we get the correct behavior, and it was missed. Will push the correct files shortly
Also, since this is supposed to be a local only potential, the cutoff is irrelevant in that case. The local potential is always splined on the original grid defined in the xml file, which always from from 0 to 10 in ppconvert I believe