Dear PySCF community

I was testing the Restricted-DHF:
for the Helium atom, with a really small basis, the solution converges:

mol = gto.M(atom='He', basis='sto-3g')
mf = scf.dhf.RDHF(mol)
mf.kernel()
converged SCF energy = -2.80786900919401

But it only works for that system with that specific basis. It does not work for the same system with a larger basis, or for other systems. For example, with water:

mol = gto.M(atom='H 0 0 0; H 0 0 2; O 0 0 1', basis='sto-3g')
mf = scf.dhf.RDHF(mol)
mf.kernel()

**/home/daniel/zquatev/zquatev/zquatev.py:50: RuntimeWarning: invalid value encountered in divide
  vcanon = v_valid / numpy.sqrt(e_valid)**

    return zquatev.solve_KR_FCSCE(self.mol, h, s)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/daniel/zquatev/zquatev/zquatev.py", line 148, in solve_KR_FCSCE
    e,v = geigh(tfock,tova,debug)
          ^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/daniel/zquatev/zquatev/zquatev.py", line 57, in geigh
    e_fock, v_fock = eigh(tfockmo, iop=1)
                     ^^^^^^^^^^^^^^^^^^^^
  File "/home/daniel/zquatev/zquatev/zquatev.py", line 9, in eigh
    err = libzquatev.eigh(fmat,eigs)
          ^^^^^^^^^^^^^^^^^^^^^^^^^^
**RuntimeError: zquatev failed**

But, keeping only with the positives eigenvalues in the geigh function, the solutions doesn't converge or the diis eigenvalues became singular:

    thresh = 1e-8
    mask = e > thresh
    e_valid = e[mask]
    v_valid = v[:, mask]

WARN:  diis singular, eigh(h) [-2.77176977e-07 -9.84229062e-12  9.13970225e-09  6.00182661e-08
  2.77713032e-07  2.27467808e-06  1.11325477e-04  4.67261717e-03
  2.90221348e+07]
SCF not converged.
SCF energy = -4.56018508433646

My environment info:

Environment Information

• OS: Linux 7.0.0-15-generic x86_64
• Python: 3.12.13
• PySCF: 2.13.0
• NumPy: 2.4.4
• SciPy: 1.17.1

Does anyone could use safely the RDHF class? If so, how have you handle this issues?

Thanks in advance!
Best

Daniel

Originally posted by @DanielBajac in #1941