RHF with H2 errors

[mtfishman]

Hi @gustavojra,

I'm finding the RHF on H2 fails:

using Fermi

@molecule {
  H  0.0 0.0 0.0
  H  0.76 0.0 0.0
}
@set basis sto-3g
wfn = @energy rhf

with the output:

================================================================================
|                                 Hartree-Fock                                 |
|                                  Module  by                                  |
|                         G.J.R Aroeira and M.M. Davis                         |
================================================================================
Collecting necessary integrals...
Done in    0.57259 s
Using GWH Guess
Molecule:

H    0.000000000000    0.000000000000    0.000000000000
H    0.760000000000    0.000000000000    0.000000000000


Charge: 0   Multiplicity: 1   
Nuclear repulsion:    0.6962858038
 Number of AOs:                            2
 Number of Doubly Occupied Orbitals:       1
 Number of Virtual Spatial Orbitals:       1
nt: 6
 Guess Energy    -1.81166646326759

 Iter.            E[RHF]         ΔE       Dᵣₘₛ        t     DIIS     damp
--------------------------------------------------------------------------------
    1     -1.1153806594  -1.115e+00   0.000e+00     0.76    false      NaN
ERROR: LoadError: ArgumentError: matrix contains Infs or NaNs
Stacktrace:
  [1] chkuplofinite(A::Matrix{Float64}, uplo::Char)
    @ LinearAlgebra.LAPACK /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/lapack.jl:109
  [2] syevr!(jobz::Char, range::Char, uplo::Char, A::Matrix{Float64}, vl::Float64, vu::Float64, il::Int64, iu::Int64, abstol::Float64)
    @ LinearAlgebra.LAPACK /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/lapack.jl:5089
  [3] #eigen!#185
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/symmetriceigen.jl:4 [inlined]
  [4] #eigen#186
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/symmetriceigen.jl:9 [inlined]
  [5] #diagonalize#20
    @ ~/.julia/packages/Fermi/phTqJ/src/Backend/Arrays.jl:189 [inlined]
  [6] macro expansion
    @ ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:138 [inlined]
  [7] macro expansion
    @ ./timing.jl:299 [inlined]
  [8] macro expansion
    @ ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:133 [inlined]
  [9] macro expansion
    @ ./timing.jl:299 [inlined]
 [10] Fermi.HartreeFock.RHF(ints::Fermi.Integrals.IntegralHelper{Float64, Fermi.Integrals.SparseERI, Fermi.Orbitals.AtomicOrbitals}, C::FermiMDArray{Float64, 2}, Λ::FermiMDArray{Float64, 2}, Alg::Fermi.HartreeFock.RHFa)
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:132
 [11] Fermi.HartreeFock.RHF(ints::Fermi.Integrals.IntegralHelper{Float64, Fermi.Integrals.SparseERI, Fermi.Orbitals.AtomicOrbitals}, Alg::Fermi.HartreeFock.RHFa)
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:32
 [12] Fermi.HartreeFock.RHF(Alg::Fermi.HartreeFock.RHFa)
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:5
 [13] Fermi.HartreeFock.RHF()
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHF.jl:81
 [14] top-level scope
    @ ~/Dropbox (Simons Foundation)/workdir/Fermi.jl/hartree_fock_H2.jl:8
 [15] include(fname::String)
    @ Base.MainInclude ./client.jl:451
 [16] top-level scope
    @ REPL[1]:1
in expression starting at /home/mfishman/Dropbox (Simons Foundation)/workdir/Fermi.jl/hartree_fock_H2.jl:8

My setup is:

(@v1.7) pkg> st Fermi
      Status `~/.julia/environments/v1.7/Project.toml`
  [9237668d] Fermi v0.4.0

julia> versioninfo()
Julia Version 1.7.1
Commit ac5cc99908 (2021-12-22 19:35 UTC)
Platform Info:
  OS: Linux (x86_64-pc-linux-gnu)
  CPU: Intel(R) Xeon(R) E-2176M  CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-12.0.1 (ORCJIT, skylake)
Environment:
  JULIA_EDITOR = vim

Maybe I'm running it wrong.

Thanks,
Matt

[gustavojra]

Hi Matt,

It is a problem with DIIS not dealing well with small systems. It does need fixing, but I have been procrastinating on it. You can simply disable DIIS if the system is small (~5 electrons or so)

@set diis false

[mtfishman]

Ok, thanks. Glad to know it's not just me doing something silly. Unfortunately that doesn't seem to help. When I run:

using Fermi

@molecule {
  H  0.0 0.0 0.0
  H  0.76 0.0 0.0
}
@set basis sto-3g
@set diis false
wfn = @energy rhf

I get what looks like the same error.

[gustavojra]

Oh, interesting. If you use a larger basis such as cc-pvdz does the error persist?

[mtfishman]

Looks like cc-pvdz works with our without DIIS.

[gustavojra]

Can you try @set oda false for the sto-3g case? That damp NaN on the first iteration does not appear for me.

[mtfishman]

Ok, that's seems to work:

using Fermi

@molecule {
  H  0.0 0.0 0.0
  H  0.76 0.0 0.0
}
@set basis sto-3g
@set diis false
@set oda false
wfn = @energy rhf

outputs:

================================================================================
|                                 Hartree-Fock                                 |
|                                  Module  by                                  |
|                         G.J.R Aroeira and M.M. Davis                         |
================================================================================
Collecting necessary integrals...
Done in    0.00009 s
Using GWH Guess
Molecule:

H    0.000000000000    0.000000000000    0.000000000000
H    0.760000000000    0.000000000000    0.000000000000


Charge: 0   Multiplicity: 1   
Nuclear repulsion:    0.6962858038
 Number of AOs:                            2
 Number of Doubly Occupied Orbitals:       1
 Number of Virtual Spatial Orbitals:       1
 Guess Energy    -1.81166646326759

 Iter.            E[RHF]         ΔE       Dᵣₘₛ        t     DIIS     damp
--------------------------------------------------------------------------------
    1     -1.1153806594  -1.115e+00   0.000e+00     0.00    false     0.00
    2     -1.1153806594   0.000e+00   0.000e+00     0.00    false     0.00
    3     -1.1153806594   0.000e+00   0.000e+00     0.00    false     0.00
    4     -1.1153806594   0.000e+00   0.000e+00     0.00    false     0.00
    5     -1.1153806594   0.000e+00   0.000e+00     0.00    false     0.00
--------------------------------------------------------------------------------
    RHF done in  0.00s
    @Final RHF Energy          -1.115380659448 Eₕ

   • Orbitals Summary

    Orbital            Energy    Occupancy
          1     -0.5703627667       ↿⇂
          2      0.6508750828         

   ✔  SCF Equations converged 😄
--------------------------------------------------------------------------------
 ⇒ Fermi Restricted Hartree--Fock Wave function
 ⋅ Basis:                  sto-3g
 ⋅ Energy:                 -1.1153806594478526
 ⋅ Occ. Spatial Orbitals:  1
 ⋅ Vir. Spatial Orbitals:  1
Convergence: ΔE => 0.00e+00 Dᵣₘₛ => 0.00e+00

With diis true it errors with:

================================================================================
|                                 Hartree-Fock                                 |
|                                  Module  by                                  |
|                         G.J.R Aroeira and M.M. Davis                         |
================================================================================
Collecting necessary integrals...
Done in    0.00009 s
Using GWH Guess
Molecule:

H    0.000000000000    0.000000000000    0.000000000000
H    0.760000000000    0.000000000000    0.000000000000


Charge: 0   Multiplicity: 1   
Nuclear repulsion:    0.6962858038
 Number of AOs:                            2
 Number of Doubly Occupied Orbitals:       1
 Number of Virtual Spatial Orbitals:       1
 Guess Energy    -1.81166646326759

 Iter.            E[RHF]         ΔE       Dᵣₘₛ        t     DIIS     damp
--------------------------------------------------------------------------------
    1     -1.1153806594  -1.115e+00   0.000e+00     0.00    false     0.00
    2     -1.1153806594   0.000e+00   0.000e+00     0.00    false     0.00
    3     -1.1153806594   0.000e+00   0.000e+00     0.00    false     0.00
ERROR: LoadError: LinearAlgebra.SingularException(2)
Stacktrace:
  [1] checknonsingular
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/factorization.jl:19 [inlined]
  [2] checknonsingular
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/factorization.jl:21 [inlined]
  [3] #lu!#146
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/lu.jl:82 [inlined]
  [4] #lu#153
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/lu.jl:279 [inlined]
  [5] lu (repeats 2 times)
    @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/lu.jl:278 [inlined]
  [6] \(A::Matrix{Float64}, B::Vector{Float64})
    @ LinearAlgebra /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/LinearAlgebra/src/generic.jl:1142
  [7] extrapolate(M::Fermi.DIIS.DIISManager{Float64, Float64}; add_res::Bool)
    @ Fermi.DIIS ~/.julia/packages/Fermi/phTqJ/src/Core/DIIS.jl:83
  [8] extrapolate
    @ ~/.julia/packages/Fermi/phTqJ/src/Core/DIIS.jl:61 [inlined]
  [9] macro expansion
    @ ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:182 [inlined]
 [10] macro expansion
    @ ./timing.jl:299 [inlined]
 [11] macro expansion
    @ ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:133 [inlined]
 [12] macro expansion
    @ ./timing.jl:299 [inlined]
 [13] Fermi.HartreeFock.RHF(ints::Fermi.Integrals.IntegralHelper{Float64, Fermi.Integrals.SparseERI, Fermi.Orbitals.AtomicOrbitals}, C::FermiMDArray{Float64, 2}, Λ::FermiMDArray{Float64, 2}, Alg::Fermi.HartreeFock.RHFa)
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:132
 [14] Fermi.HartreeFock.RHF(ints::Fermi.Integrals.IntegralHelper{Float64, Fermi.Integrals.SparseERI, Fermi.Orbitals.AtomicOrbitals}, Alg::Fermi.HartreeFock.RHFa)
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:32
 [15] Fermi.HartreeFock.RHF(Alg::Fermi.HartreeFock.RHFa)
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHFa.jl:5
 [16] Fermi.HartreeFock.RHF()
    @ Fermi.HartreeFock ~/.julia/packages/Fermi/phTqJ/src/Methods/HartreeFock/RHF/RHF.jl:81
 [17] top-level scope
    @ ~/Dropbox (Simons Foundation)/workdir/Fermi.jl/hartree_fock_H2.jl:10
 [18] include(fname::String)
    @ Base.MainInclude ./client.jl:451
 [19] top-level scope
    @ REPL[6]:1
in expression starting at /home/mfishman/Dropbox (Simons Foundation)/workdir/Fermi.jl/hartree_fock_H2.jl:10

[gustavojra]

Okay, great. It's a good reminder to add numerical checks for these extrapolations schemes. They should be well behaved for larger systems though.