Potential problem on exponent of carbon A5ZP basis set

[ajz34]

Hi BSE devevlopers!

I found that for carbon atom, the exponent of diffusion basis function of Jorge-A5ZP d shell is 0.007196; this exponent is probably too diffuse, comparing to 0.043962 of boron or 0.09339 of nitrogen.

https://github.com/MolSSI-BSE/basis_set_exchange/blob/f8c059b9ecf34571b358b04a5303b4b4bd8c9268/basis_set_exchange/data/jorge/A5ZP.1.json#L1933

I guess the exponent of carbon diffuse d shell is probably 0.07196.

I also confirmed that the basis set data of basis_set_exchange is the same to that obtained from https://qcgv.ufes.br/. So it's probably not a bug of basis_set_exchange, but potential typo of A5ZP basis set itself.

[susilehtola]

It might not even be a typo. Did you contact Jorge?

[ajz34]

@susilehtola Currently not contacted. I'll mail to Jorge on this topic.

[ajz34]

@susilehtola After contacting to Jorge, I think that it is actually not a typo. You are correct and I was less considerated on this topic 😿

I can verify that exponents of diffuse GTO primitives of (Jorge-)ADZP (d shell) and ATZP (d, f shells) for carbon could be reproduced by minimizing ROMP2 atomic anion energy (described by 10.1016/j.theochem.2004.11.037). Though I haven't fully reproduced exponents of AQZP and A5ZP by optimization, d shell exponent 0.007196 of carbon gives lower ROMP2 energy (-37.8191569 a.u. with 2 frozen electrons) than 0.07196 (-37.8146132 a.u.). So it's reasonable to apply 0.007196 as exponent for d symmetry as diffusion primitive of carbon atom, from this perspective.

[susilehtola]

When I repeat the ROMP2 calculations with Gaussian, I get -37.8191670 for the original exponent and -37.8146132 for the exponent multiplied by 10. The latter number matches Jorge's, but there's a discrepancy in the first one, which is already indicative of a problem. Namely, examining the log file it becomes clear that the ROHF calculation didn't even converge but Gaussian went ahead and computed a UMP2 energy. It looks like the exponent comes out so small because of this convergence failure, which makes MP2 completely unreliable.

Repeating the calculation with UMP2, the optimal exponent turns out to come out around 0.06, as one would expect based on the boron and nitrogen calculation.

Similar issues likely affect all of Jorge's basis sets optimized with ROMP2 in Gaussian.

[susilehtola]

Here's the scan for proof

[susilehtola]

I have emailed Jorge about this issue.

[susilehtola]

I also emailed Gaussian tech support. The answer came pretty quickly: the DIIS gets confused by small exponents, but the calculations run correctly if you set Use=L506 in the route section as is also directed in https://gaussian.com/scf/