Very different H3 molecule ROHF energy with near identical geometry

Hello developers,

I have been using OpenFermion (OF), which is a package for quantum computing. OF calls PSI4 for the ab initio calculation and parse the output to give a second-quantized Hamiltonian.

OF failed in many cases for the H3 neutral molecule ground state energy, in the sense that when you dissociate H3 as an equilateral triangle gradually, the HF energy would change abruptly after certain dissociation distance. I suspect it likely has something to do with the fact that PSI4 gives very unstable ground state energy at these distances. Below is the minimal input to reproduce the error.

I will show input 1, input 2, and then their outputs. I truncated output due to the length. But they should be easily reproduced using psi4 input.
The observation is that, they have very different ROHF energy but the same FCI energy.
I don’t know how to explain that, and I don’t know which ROHF number to trust.
I want the energy to be consistent (meaning roughly the same) for H3 with near-identical geometry so that the H3 equilateral dissociation curve is smooth. And I hope there is a way to obtain the correct ROHF energy reliably.

Input 1

molecule mol {
    H 0.0        3.1000000000000014 0.0
    H 0.0 0.0 0.0
    H        2.6846787517317607        1.5500000000000012 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
    fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf false
    scf_type pk
    maxiter 100
    num_amps_print 1e6
    r_convergence 1e-6
    d_convergence 1e-6
    e_convergence 1e-6
    ints_tolerance EQUALITY_TOLERANCE
    damping_percentage 0
}
set qc_module detci

hf_energy= energy('scf')
fci_energy = energy('fci')

Input 2


molecule mol {
    H 0.0        3.1000000000000001 0.0
    H 0.0 0.0 0.0
    H        2.6846787517317598        1.5500000000000005 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
    fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf false
    scf_type pk
    maxiter 100
    num_amps_print 1e6
    r_convergence 1e-6
    d_convergence 1e-6
    e_convergence 1e-6
    ints_tolerance EQUALITY_TOLERANCE
    damping_percentage 0
}
set qc_module detci

hf_energy= energy('scf')
fci_energy = energy('fci')

Output 1


    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.4a2.dev1055 

                         Git: Rev {loriab-patch-1-hokru-soln} 921c1d2 



         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H           -0.000000000000     0.000000000000     1.789785834488     1.007825032230
         H           -1.550000000000     0.000000000000    -0.894892917244     1.007825032230
         H            1.550000000000     0.000000000000    -0.894892917244     1.007825032230

  Running in c1 symmetry.

  ==> Pre-Iterations <==

  SCF Guess: Superposition of Atomic Densities via on-the-fly atomic UHF (no occupation information).

   -------------------------
    Irrep   Nso     Nmo    
   -------------------------
     A          3       3 
   -------------------------
    Total       3       3
   -------------------------

  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter SAD:    -0.81895913300920   -8.18959e-01   0.00000e+00 
   @ROHF iter   1:    -1.02462826495174   -2.05669e-01   4.67186e-02 DIIS
   @ROHF iter   2:    -0.81904807260243    2.05580e-01   1.78670e-02 DIIS
   @ROHF iter   3:    -0.82429905280403   -5.25098e-03   2.72163e-02 DIIS
   @ROHF iter   4:    -0.81733207095513    6.96698e-03   2.29653e-02 DIIS
   @ROHF iter   5:    -0.81900780327629   -1.67573e-03   1.79226e-02 DIIS
   @ROHF iter   6:    -0.81948849457445   -4.80691e-04   1.81392e-02 DIIS
   @ROHF iter   7:    -0.81825999019787    1.22850e-03   1.77550e-02 DIIS
   @ROHF iter   8:    -0.82558019562478   -7.32021e-03   2.34476e-02 DIIS
   @ROHF iter   9:    -0.80773370032587    1.78465e-02   1.34632e-02 DIIS
   @ROHF iter  10:    -0.81507682881535   -7.34313e-03   1.50698e-02 DIIS
   @ROHF iter  11:    -1.09347817715635   -2.78401e-01   2.70754e-02 DIIS
   @ROHF iter  12:    -0.80333433494676    2.90144e-01   8.18897e-03 DIIS
   @ROHF iter  13:    -0.88110669412968   -7.77724e-02   4.43919e-02 DIIS
   @ROHF iter  14:    -1.10680822410247   -2.25702e-01   1.68274e-02 DIIS
   @ROHF iter  15:    -1.10089321361413    5.91501e-03   2.13647e-02 DIIS
   @ROHF iter  16:    -0.87413876923508    2.26754e-01   4.22718e-02 DIIS
   @ROHF iter  17:    -0.80291475256590    7.12240e-02   6.42437e-03 DIIS
   @ROHF iter  18:    -0.82040461366209   -1.74899e-02   2.12211e-02 DIIS
   @ROHF iter  19:    -0.81069052393353    9.71409e-03   1.16020e-02 DIIS
   @ROHF iter  20:    -0.80978189471290    9.08629e-04   1.01562e-02 DIIS
   @ROHF iter  21:    -0.80863392256624    1.14797e-03   8.51342e-03 DIIS
   @ROHF iter  22:    -0.80468403971157    3.94988e-03   1.07202e-02 DIIS
   @ROHF iter  23:    -1.10131664166576   -2.96633e-01   2.19974e-02 DIIS
   @ROHF iter  24:    -0.81146017154109    2.89856e-01   1.18683e-02 DIIS
   @ROHF iter  25:    -0.81284187287457   -1.38170e-03   1.33734e-02 DIIS
   @ROHF iter  26:    -1.08474128934510   -2.71899e-01   3.09169e-02 DIIS
   @ROHF iter  27:    -1.10658531285378   -2.18440e-02   1.70679e-02 DIIS
   @ROHF iter  28:    -1.10614793125730    4.37382e-04   1.77240e-02 DIIS
   @ROHF iter  29:    -1.10236290839220    3.78502e-03   2.07016e-02 DIIS
   @ROHF iter  30:    -1.10618344981553   -3.82054e-03   1.72371e-02 DIIS
   @ROHF iter  31:    -0.80394958725847    3.02234e-01   9.27080e-03 DIIS
   @ROHF iter  32:    -0.80291228980406    1.03730e-03   8.42285e-03 DIIS
   @ROHF iter  33:    -0.80875165935548   -5.83937e-03   1.60142e-02 DIIS
   @ROHF iter  34:    -0.80840037415664    3.51285e-04   1.55749e-02 DIIS
   @ROHF iter  35:    -0.80888287423803   -4.82500e-04   1.60143e-02 DIIS
   @ROHF iter  36:    -0.80819081935446    6.92055e-04   1.52825e-02 DIIS
   @ROHF iter  37:    -0.80856912378300   -3.78304e-04   1.56786e-02 DIIS
   @ROHF iter  38:    -0.80830050475233    2.68619e-04   1.53962e-02 DIIS
   @ROHF iter  39:    -0.80847243960933   -1.71935e-04   1.55771e-02 DIIS
   @ROHF iter  40:    -0.80835710334830    1.15336e-04   1.54559e-02 DIIS
   @ROHF iter  41:    -0.80843254264875   -7.54393e-05   1.55352e-02 DIIS
   @ROHF iter  42:    -0.80838245862410    5.00840e-05   1.54826e-02 DIIS
   @ROHF iter  43:    -0.80841540071944   -3.29421e-05   1.55172e-02 DIIS
   @ROHF iter  44:    -0.80839360255846    2.17982e-05   1.54943e-02 DIIS
   @ROHF iter  45:    -0.80840797016406   -1.43676e-05   1.55094e-02 DIIS
   @ROHF iter  46:    -0.80839847568326    9.49448e-06   1.54994e-02 DIIS
   @ROHF iter  47:    -0.80840473922804   -6.26354e-06   1.55060e-02 DIIS
   @ROHF iter  48:    -0.80840060250417    4.13672e-06   1.55016e-02 DIIS
   @ROHF iter  49:    -0.80840333256071   -2.73006e-06   1.55045e-02 DIIS
   @ROHF iter  50:    -0.80840152996258    1.80260e-06   1.55026e-02 DIIS
   @ROHF iter  51:    -0.80840271979607   -1.18983e-06   1.55039e-02 DIIS
   @ROHF iter  52:    -0.80840193426060    7.85535e-07   1.55030e-02 DIIS
   @ROHF iter  53:    -0.80840245280273   -5.18542e-07   1.55036e-02 DIIS
   @ROHF iter  54:    -0.80840211047506    3.42328e-07   1.55032e-02 DIIS
   @ROHF iter  55:    -0.80840233645692   -2.25982e-07   1.55035e-02 DIIS
   @ROHF iter  56:    -0.80840218727221    1.49185e-07   1.55033e-02 DIIS
   @ROHF iter  57:    -0.80840228575545   -9.84832e-08   1.55034e-02 DIIS
   @ROHF iter  58:    -0.80840222074147    6.50140e-08   1.55033e-02 DIIS
   @ROHF iter  59:    -0.80840226366015   -4.29187e-08   1.55034e-02 DIIS
   @ROHF iter  60:    -0.80840223532716    2.83330e-08   1.55034e-02 DIIS
   @ROHF iter  61:    -0.80840225403157   -1.87044e-08   1.55034e-02 DIIS
   @ROHF iter  62:    -0.80840224168383    1.23477e-08   1.55034e-02 DIIS
   @ROHF iter  63:    -0.80840224983506   -8.15123e-09   1.55034e-02 DIIS
   @ROHF iter  64:    -0.80840224445432    5.38074e-09   1.55034e-02 DIIS
   @ROHF iter  65:    -0.80840224800552   -3.55120e-09   1.55034e-02 DIIS
   @ROHF iter  66:    -0.80840224566177    2.34375e-09   1.55034e-02 DIIS
   @ROHF iter  67:    -0.80840224720874   -1.54696e-09   1.55034e-02 DIIS
   @ROHF iter  68:    -0.80840224618751    1.02123e-09   1.55034e-02 DIIS
   @ROHF iter  69:    -0.80840224686208   -6.74571e-10   1.55034e-02 DIIS
   @ROHF iter  70:    -0.80840224641602    4.46056e-10   1.55034e-02 DIIS
   @ROHF iter  71:    -0.80840224671100   -2.94981e-10   1.55034e-02 DIIS
   @ROHF iter  72:    -0.80840224651689    1.94111e-10   1.55034e-02 DIIS
   @ROHF iter  73:    -0.80840224664395   -1.27062e-10   1.55034e-02 DIIS
   @ROHF iter  74:    -0.80840224656024    8.37121e-11   1.55034e-02 DIIS
   @ROHF iter  75:    -0.80840224661622   -5.59777e-11   1.55034e-02 DIIS
   @ROHF iter  76:    -0.80840224657924    3.69769e-11   1.55034e-02 DIIS
   @ROHF iter  77:    -0.80840224660323   -2.39884e-11   1.55034e-02 DIIS
   @ROHF iter  78:    -0.80840224658712    1.61122e-11   1.55034e-02 DIIS
   @ROHF iter  79:    -0.80840224659835   -1.12268e-11   1.55034e-02 DIIS
   @ROHF iter  80:    -0.80840224659048    7.86970e-12   1.55034e-02 DIIS
   @ROHF iter  81:    -0.80840224659568   -5.20939e-12   1.55034e-02 DIIS
   @ROHF iter  82:    -0.80840224659236    3.32090e-12   1.55034e-02 DIIS
   @ROHF iter  83:    -0.80840224659425   -1.88227e-12   1.55034e-02 DIIS
   @ROHF iter  84:    -0.80840224659324    1.00431e-12   1.55034e-02 DIIS
   @ROHF iter  85:    -0.80840224659458   -1.33982e-12   1.55034e-02 DIIS
   @ROHF iter  86:    -0.80840224659274    1.84541e-12   1.55034e-02 DIIS
   @ROHF iter  87:    -0.80840224659452   -1.78635e-12   1.55034e-02 DIIS
   @ROHF iter  88:    -0.80840224659361    9.13714e-13   1.55034e-02 DIIS
   @ROHF iter  89:    -0.80840224659377   -1.59872e-13   1.55034e-02 DIIS
   @ROHF iter  90:    -0.80840224659331    4.61853e-13   1.55034e-02 DIIS
   @ROHF iter  91:    -0.80840224659378   -4.69846e-13   1.55034e-02 DIIS
   @ROHF iter  92:    -0.80840224659415   -3.72813e-13   1.55034e-02 DIIS
   @ROHF iter  93:    -0.80840224659355    6.03739e-13   1.55034e-02 DIIS
   @ROHF iter  94:    -0.80840224659341    1.37668e-13   1.55034e-02 DIIS
   @ROHF iter  95:    -0.80840224659392   -5.12479e-13   1.55034e-02 DIIS
   @ROHF iter  96:    -0.80840224659348    4.36096e-13   1.55034e-02 DIIS
   @ROHF iter  97:    -0.80840224659358   -9.90319e-14   1.55034e-02 DIIS
   @ROHF iter  98:    -0.80840224659367   -8.57092e-14   1.55034e-02 DIIS
   @ROHF iter  99:    -0.80840224659347    1.95843e-13   1.55034e-02 DIIS
   @ROHF iter 100:    -0.80840224659401   -5.37570e-13   1.55034e-02 DIIS

PsiException: Could not converge SCF iterations in 100 iterations.

  Energy and/or wave function did not converge, but proceeding anyway.


  ==> Post-Iterations <==

    Orbital Energies [Eh]
    ---------------------

    Doubly Occupied:                                                      

       1A     -0.207067  

    Singly Occupied:                                                      

       2A     -0.042712  

    Virtual:                                                              

       3A      0.011580  

    Final Occupation by Irrep:
              A 
    DOCC [     1 ]
    SOCC [     1 ]

  @ROHF Final Energy:    -0.80840224659401

   => Energetics <=

    Nuclear Repulsion Energy =              0.5121069780677417
    One-Electron Energy =                  -2.4329519292604904
    Two-Electron Energy =                   1.1124427045987375
    Total Energy =                         -0.8084022465940113

 



         ---------------------------------------------------------
                          Configuration Interaction
                            (a 'D E T C I' module)

                 C. David Sherrill, Daniel G. A. Smith, and
                              Matt L. Leininger
         ---------------------------------------------------------


    H0 Block Eigenvalue =  -1.39998699

    Simultaneous Expansion Method (Block Davidson Method)
    Using 1 initial trial vectors

     Iter   Root       Total Energy       Delta E      C RMS

   @CI  0:     0     -1.399986986556   -1.9121E+00   1.4676E-15  
    Warning: Norm of correction (root 0) is < 1.0E-13
   @CI  1:     0     -1.399986986556   0.0000E+00   1.3714E-15 c

   ==> Energetics <==

    SCF energy =           -0.808402246594011
    Total CI energy =      -1.399986986555535

   ==> FCI root 0 information <==

    FCI Root 0 energy =    -1.399986986555535

   The 9 most important determinants:

    *   1   -0.619973  (    2,    0)  1AB 2AA 3AA 
    *   2    0.459916  (    0,    1)  1AA 2AX 
    *   3   -0.435287  (    1,    2)  1AA 3AX 
    *   4   -0.382980  (    1,    1)  1AA 2AB 3AA 
    *   5    0.236993  (    0,    2)  1AA 2AA 3AB 
    *   6   -0.094055  (    1,    0)  1AX 3AA 
    *   7    0.051133  (    2,    1)  2AX 3AA 
    *   8    0.018365  (    0,    0)  1AX 2AA
		

Output 2

    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.4a2.dev1055 

                         Git: Rev {loriab-patch-1-hokru-soln} 921c1d2 






         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H            0.000000000000     0.000000000000     1.789785834488     1.007825032230
         H           -1.550000000000     0.000000000000    -0.894892917244     1.007825032230
         H            1.550000000000     0.000000000000    -0.894892917244     1.007825032230

  Running in c1 symmetry.


  

 


  ==> Pre-Iterations <==

  SCF Guess: Superposition of Atomic Densities via on-the-fly atomic UHF (no occupation information).

   -------------------------
    Irrep   Nso     Nmo    
   -------------------------
     A          3       3 
   -------------------------
    Total       3       3
   -------------------------

  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter SAD:    -0.81895913300919   -8.18959e-01   0.00000e+00 
   @ROHF iter   1:    -1.02462826495174   -2.05669e-01   2.47354e-02 DIIS
   @ROHF iter   2:    -0.95016024840700    7.44680e-02   5.02761e-02 DIIS
   @ROHF iter   3:    -0.99605770375080   -4.58975e-02   4.35380e-02 DIIS
   @ROHF iter   4:    -1.04967524373170   -5.36175e-02   8.35149e-03 DIIS
   @ROHF iter   5:    -1.05474365190057   -5.06841e-03   1.85273e-02 DIIS
   @ROHF iter   6:    -1.04863546782488    6.10818e-03   8.70364e-04 DIIS
   @ROHF iter   7:    -1.04865114053343   -1.56727e-05   4.15913e-04 DIIS
   @ROHF iter   8:    -1.04864685121106    4.28932e-06   1.92684e-04 DIIS
   @ROHF iter   9:    -1.04864709143170   -2.40221e-07   2.37095e-05 DIIS
   @ROHF iter  10:    -1.04864710193089   -1.04992e-08   1.38084e-07 DIIS
  Energy and wave function converged.


  ==> Post-Iterations <==

    Orbital Energies [Eh]
    ---------------------

    Doubly Occupied:                                                      

       1A     -0.185028  

    Singly Occupied:                                                      

       2A     -0.068033  

    Virtual:                                                              

       3A      0.011745  

    Final Occupation by Irrep:
              A 
    DOCC [     1 ]
    SOCC [     1 ]

  @ROHF Final Energy:    -1.04864710193089

   => Energetics <=

    Nuclear Repulsion Energy =              0.5121069780677421
    One-Electron Energy =                  -2.4452074588905082
    Two-Electron Energy =                   0.8844533788918788
    Total Energy =                         -1.0486471019308872

Computation Completed


Properties will be evaluated at   0.000000,   0.000000,   0.000000 [a0]

Properties computed using the SCF density matrix

  Nuclear Dipole Moment: [e a0]
     X:     0.0000      Y:     0.0000      Z:     0.0000

  Electronic Dipole Moment: [e a0]
     X:     0.0000      Y:     0.0000      Z:     0.0809

  Dipole Moment: [e a0]
     X:     0.0000      Y:     0.0000      Z:     0.0809     Total:     0.0809

  Dipole Moment: [D]
     X:     0.0000      Y:     0.0000      Z:     0.2057     Total:     0.2057


  ==> Molecular Orbitals <==

                           1            2            3

 1    H1 s0         0.6941755   -0.0000000    0.7201039
 2    H2 s0         0.4952751    0.7123810   -0.4976394
 3    H3 s0         0.4952751   -0.7123810   -0.4976394

            Ene    -0.1850277   -0.0680326    0.0117453
            Sym             A            A            A
            Occ             2            1            0


*** tstop() called on wireless-10-104-132-196.umd.edu at Mon Apr 18 12:20:00 2022
Module time:
	user time   =       0.22 seconds =       0.00 minutes
	system time =       0.01 seconds =       0.00 minutes
	total time  =          1 seconds =       0.02 minutes
Total time:
	user time   =       0.45 seconds =       0.01 minutes
	system time =       0.03 seconds =       0.00 minutes
	total time  =          1 seconds =       0.02 minutes
 MINTS: Wrapper to libmints.
   by Justin Turney

   Calculation information:
      Number of threads:                 1
      Number of atoms:                   3
      Number of AO shells:               3
      Number of SO shells:               3
      Number of primitives:              9
      Number of atomic orbitals:         3
      Number of basis functions:         3

      Number of irreps:                  1
      Integral cutoff                 0.00e+00
      Number of functions per irrep: [   3 ]

 OEINTS: Overlap, kinetic, potential, dipole, and quadrupole integrals
         stored in file 35.

      Computing two-electron integrals...done
      Computed 21 non-zero two-electron integrals.
        Stored in file 33.


         ---------------------------------------------------------
                          Configuration Interaction
                            (a 'D E T C I' module)

                 C. David Sherrill, Daniel G. A. Smith, and
                              Matt L. Leininger
         ---------------------------------------------------------

Warning: iopen=0,opentype!=closed. Making iopen=1

   

   ==> CI Orbital and Space information <==

   ------------------------------------
               Space    Total     A
   ------------------------------------
                 Nso        3     3
                 Nmo        3     3
               Ndocc        1     1
               Nsocc        1     1
   ------------------------------------
              CI Spaces
   ------------------------------------
        Dropped DOCC        0     0
              Active        3     3
        Dropped UOCC        0     0
   ------------------------------------

   ==> Setting up CI strings <==

    There are 3 alpha and 3 beta strings
    The CI space requires 9 (9.00E+00) determinants and 1 blocks

   ==> Transforming CI integrals <==

	Presorting SO-basis two-electron integrals.
	Sorting File: SO Ints (nn|nn) nbuckets = 1
	Constructing frozen core operators
	Starting first half-transformation.
	Sorting half-transformed integrals.
	First half integral transformation complete.
	Starting second half-transformation.
	Two-electron integral transformation complete.

   ==> Starting CI iterations <==

    H0 Block Eigenvalue =  -1.39998699

    Simultaneous Expansion Method (Block Davidson Method)
    Using 1 initial trial vectors

     Iter   Root       Total Energy       Delta E      C RMS

   @CI  0:     0     -1.399986986556   -1.9121E+00   1.3629E-15  
    Warning: Norm of correction (root 0) is < 1.0E-13
   @CI  1:     0     -1.399986986556   0.0000E+00   2.7565E-15 c

   ==> Energetics <==

    SCF energy =           -1.048647101930887
    Total CI energy =      -1.399986986555534

   ==> FCI root 0 information <==

    FCI Root 0 energy =    -1.399986986555534

   The 9 most important determinants:

    *   1    0.512161  (    1,    0)  1AX 3AA 
    *   2   -0.482488  (    1,    2)  1AA 3AX 
    *   3    0.474268  (    0,    1)  1AA 2AX 
    *   4   -0.461834  (    2,    1)  2AX 3AA 
    *   5    0.164880  (    0,    0)  1AX 2AA 
    *   6    0.151250  (    2,    2)  2AA 3AX 
    *   7    0.105207  (    1,    1)  1AA 2AB 3AA 
    *   8    0.055040  (    0,    2)  1AA 2AA 3AB 
    *   9    0.050167  (    2,    0)  1AB 2AA 3AA 

		

Your copy of Psi4 isn’t able to solve the ROHF equations for the first input. Normally, Psi4 would raise an error telling you this, but you specifically set fail_on_maxiter False, which tells Psi4 not to raise such an error. This is what

PsiException: Could not converge SCF iterations in 100 iterations.

  Energy and/or wave function did not converge, but proceeding anyway.

means.

The question then is “how do I get Psi4 to solve the ROHF equations correctly?” I can’t reproduce your ROHF problems on my copy of Psi. What is Git: Rev {loriab-patch-1-hokru-soln} 921c1d2 appearing in your Psi version information? This looks like you have made custom edits to your version of Psi. (Which is, by the way, outdated. Please upgrade to v1.5.)

Thanks for looking into this. Yes, I see that the energy different is probably due to the failed convergence. It would be great to know how to converge the ROHF calculation reliably.

I have no idea what that psi4 version is. Anyway, I tried to create a new conda environment with the psi4 1.5 version, and I still got the wrong result. See below.

Output 1


    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.5 release

                         Git: Rev {HEAD} e9f4d6d 
    Psi4 started on: Monday, 18 April 2022 11:42PM

    
  ==> Input File <==

--------------------------------------------------------------------------

molecule mol {
    H 0.0        3.1000000000000014 0.0
    H 0.0 0.0 0.0
    H        2.6846787517317607        1.5500000000000012 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
    fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf false
    scf_type pk
    maxiter 100
    num_amps_print 1e6
    r_convergence 1e-6
    d_convergence 1e-6
    e_convergence 1e-6
    ints_tolerance EQUALITY_TOLERANCE
    damping_percentage 0
}
set qc_module detci

fci_energy = energy('fci')

--------------------------------------------------------------------------


         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H           -0.000000000000     0.000000000000     1.789785834488     1.007825032230
         H           -1.550000000000     0.000000000000    -0.894892917244     1.007825032230
         H            1.550000000000     0.000000000000    -0.894892917244     1.007825032230



  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter SAD:    -0.81895913300920   -8.18959e-01   0.00000e+00 
   @ROHF iter   1:    -1.02462826495174   -2.05669e-01   3.38204e-02 DIIS
   @ROHF iter   2:    -0.84909105626904    1.75537e-01   3.33176e-02 DIIS
   @ROHF iter   3:    -0.93431924896954   -8.52282e-02   5.09467e-02 DIIS
   @ROHF iter   4:    -0.93529395637096   -9.74707e-04   5.11293e-02 DIIS
   @ROHF iter   5:    -0.97345426880694   -3.81603e-02   5.20987e-02 DIIS
   @ROHF iter   6:    -1.00440933892290   -3.09551e-02   4.98416e-02 DIIS
   @ROHF iter   7:    -1.05286751801099   -4.84582e-02   4.16778e-02 DIIS
   @ROHF iter   8:    -1.10423070053484   -5.13632e-02   1.97773e-02 DIIS
   @ROHF iter   9:    -1.08909353775849    1.51372e-02   2.71352e-02 DIIS
   @ROHF iter  10:    -1.09421195679198   -5.11842e-03   2.56803e-02 DIIS
   @ROHF iter  11:    -1.09687382141191   -2.66186e-03   2.48129e-02 DIIS
   @ROHF iter  12:    -1.10710324624622   -1.02294e-02   1.72284e-02 DIIS
   @ROHF iter  13:    -1.11423214680511   -7.12890e-03   8.01588e-03 DIIS
   @ROHF iter  14:    -1.11493792690878   -7.05780e-04   6.36415e-03 DIIS
   @ROHF iter  15:    -1.11591860463621   -9.80678e-04   2.53931e-03 DIIS
   @ROHF iter  16:    -1.11610279362323   -1.84189e-04   1.23276e-04 DIIS
   @ROHF iter  17:    -1.11610323192736   -4.38304e-07   1.35890e-05 DIIS
   @ROHF iter  18:    -1.11610323717867   -5.25131e-09   1.27068e-06 DIIS
   @ROHF iter  19:    -1.11610323722485   -4.61773e-11   4.64755e-09 DIIS
  Energy and wave function converged.


  ==> Post-Iterations <==

    Orbital Energies [Eh]
    ---------------------

    Doubly Occupied:                                                      

       1A     -0.176846  

    Singly Occupied:                                                      

       2A     -0.076665  

    Virtual:                                                              

       3A      0.013252  

    Final Occupation by Irrep:
              A 
    DOCC [     1 ]
    SOCC [     1 ]

  @ROHF Final Energy:    -1.11610323722485

   => Energetics <=

    Nuclear Repulsion Energy =              0.5121069780677417
    One-Electron Energy =                  -2.4396854336095521
    Two-Electron Energy =                   0.8114752183169628
    Total Energy =                         -1.1161032372248476



         ---------------------------------------------------------
                          Configuration Interaction
                            (a 'D E T C I' module)

                 C. David Sherrill, Daniel G. A. Smith, and
                              Matt L. Leininger
         ---------------------------------------------------------

     Iter   Root       Total Energy       Delta E      C RMS

   @CI  0:     0     -1.399986986556   -1.9121E+00   2.1088E-15  
    Warning: Norm of correction (root 0) is < 1.0E-13
   @CI  1:     0     -1.399986986556   -2.2204E-16   4.3940E-15 c

   ==> Energetics <==

    SCF energy =           -1.116103237224848
    Total CI energy =      -1.399986986555537

   ==> FCI root 0 information <==

    FCI Root 0 energy =    -1.399986986555537

   The 9 most important determinants:

    *   1   -0.727652  (    0,    0)  1AX 2AA 
    *   2   -0.683226  (    2,    2)  2AA 3AX 
    *   3    0.046087  (    1,    2)  1AA 3AX 
    *   4    0.029209  (    1,    1)  1AA 2AB 3AA 
    *   5   -0.017647  (    0,    1)  1AA 2AX 
    *   6    0.014633  (    0,    2)  1AA 2AA 3AB 
    *   7    0.014576  (    2,    0)  1AB 2AA 3AA 
    *   8   -0.002984  (    1,    0)  1AX 3AA 
    *   9 

Output 2


    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.5 release

                         Git: Rev {HEAD} e9f4d6d 


  ==> Input File <==

--------------------------------------------------------------------------

molecule mol {
    H 0.0        3.1000000000000001 0.0
    H 0.0 0.0 0.0
    H        2.6846787517317598        1.5500000000000005 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
    fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf false
    scf_type pk
    maxiter 100
    num_amps_print 1e6
    r_convergence 1e-6
    d_convergence 1e-6
    e_convergence 1e-6
    ints_tolerance EQUALITY_TOLERANCE
    damping_percentage 0
}
set qc_module detci

fci_energy = energy('fci')

-------------------------------------------------------------------------


         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H            0.000000000000     0.000000000000     1.789785834488     1.007825032230
         H           -1.550000000000     0.000000000000    -0.894892917244     1.007825032230
         H            1.550000000000     0.000000000000    -0.894892917244     1.007825032230


  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter SAD:    -0.81895913300919   -8.18959e-01   0.00000e+00 
   @ROHF iter   1:    -1.02462826495174   -2.05669e-01   4.51957e-02 DIIS
   @ROHF iter   2:    -0.82090605667124    2.03722e-01   1.91612e-02 DIIS
   @ROHF iter   3:    -0.82927197929585   -8.36592e-03   2.97079e-02 DIIS
   @ROHF iter   4:    -0.82070787467750    8.56410e-03   2.50730e-02 DIIS
   @ROHF iter   5:    -0.81925047695905    1.45740e-03   1.85753e-02 DIIS
   @ROHF iter   6:    -0.82077812038270   -1.52764e-03   1.95931e-02 DIIS
   @ROHF iter   7:    -1.09284131122808   -2.72063e-01   2.53834e-02 DIIS
   @ROHF iter   8:    -1.10785015758679   -1.50088e-02   1.59352e-02 DIIS
   @ROHF iter   9:    -1.11203216304911   -4.18201e-03   1.15992e-02 DIIS
   @ROHF iter  10:    -1.09316352969779    1.88686e-02   2.52128e-02 DIIS
   @ROHF iter  11:    -1.10631442540521   -1.31509e-02   1.77531e-02 DIIS
   @ROHF iter  12:    -0.98729149417412    1.19023e-01   5.15424e-02 DIIS
   @ROHF iter  13:    -1.11126871783768   -1.23977e-01   1.28433e-02 DIIS
   @ROHF iter  14:    -0.83433214600254    2.76937e-01   3.20764e-02 DIIS
   @ROHF iter  15:    -0.80357460593966    3.07575e-02   8.55246e-03 DIIS
   @ROHF iter  16:    -1.11538193802184   -3.11807e-01   4.98225e-03 DIIS
   @ROHF iter  17:    -0.80158713171368    3.13795e-01   2.36650e-03 DIIS
   @ROHF iter  18:    -0.81062931213552   -9.04218e-03   1.09740e-02 DIIS
   @ROHF iter  19:    -0.81104948525099   -4.20173e-04   1.14253e-02 DIIS
   @ROHF iter  20:    -0.81030475612301    7.44729e-04   1.06135e-02 DIIS
   @ROHF iter  21:    -0.81039212773282   -8.73716e-05   1.07127e-02 DIIS
   @ROHF iter  22:    -0.80926447864460    1.12765e-03   9.36045e-03 DIIS
   @ROHF iter  23:    -0.80901885187203    2.45627e-04   9.66994e-03 DIIS
   @ROHF iter  24:    -0.89413078248463   -8.51119e-02   4.65685e-02 DIIS
   @ROHF iter  25:    -0.81452447920127    7.96063e-02   1.54775e-02 DIIS
   @ROHF iter  26:    -0.81348032799503    1.04415e-03   1.43331e-02 DIIS
   @ROHF iter  27:    -0.81333789978114    1.42428e-04   1.41810e-02 DIIS
   @ROHF iter  28:    -0.81331722331868    2.06765e-05   1.41592e-02 DIIS
   @ROHF iter  29:    -0.81331436102410    2.86229e-06   1.41562e-02 DIIS
   @ROHF iter  30:    -0.81331397541222    3.85612e-07   1.41558e-02 DIIS
   @ROHF iter  31:    -0.81331392401341    5.13988e-08   1.41557e-02 DIIS
   @ROHF iter  32:    -0.81331391718593    6.82748e-09   1.41557e-02 DIIS
   @ROHF iter  33:    -0.81331391628145    9.04477e-10   1.41557e-02 DIIS
   @ROHF iter  34:    -0.81331391616131    1.20139e-10   1.41557e-02 DIIS
   @ROHF iter  35:    -0.81331391614398    1.73328e-11   1.41557e-02 DIIS
   @ROHF iter  36:    -0.81331391614097    3.00826e-12   1.41557e-02 DIIS
   @ROHF iter  37:    -0.81331391614267   -1.70264e-12   1.41557e-02 DIIS
   @ROHF iter  38:    -0.81331391614287   -2.00284e-13   1.41557e-02 DIIS
   @ROHF iter  39:    -0.81331391614190    9.72555e-13   1.41557e-02 DIIS
   @ROHF iter  40:    -0.81331391614344   -1.54121e-12   1.41557e-02 DIIS
   @ROHF iter  41:    -0.81331391614399   -5.44897e-13   1.41557e-02 DIIS
   @ROHF iter  42:    -0.81331391614335    6.33493e-13   1.41557e-02 DIIS
   @ROHF iter  43:    -0.81331391614420   -8.46434e-13   1.41557e-02 DIIS
   @ROHF iter  44:    -0.81331391614152    2.67741e-12   1.41557e-02 DIIS
   @ROHF iter  45:    -0.81331391614085    6.75904e-13   1.41557e-02 DIIS
   @ROHF iter  46:    -0.81331391614290   -2.05702e-12   1.41557e-02 DIIS
   @ROHF iter  47:    -0.81331391614360   -6.98996e-13   1.41557e-02 DIIS
   @ROHF iter  48:    -0.81331391614332    2.82663e-13   1.41557e-02 DIIS
   @ROHF iter  49:    -0.81331391614364   -3.15969e-13   1.41557e-02 DIIS
   @ROHF iter  50:    -0.81331391614310    5.31575e-13   1.41557e-02 DIIS
   @ROHF iter  51:    -0.81331391614289    2.12719e-13   1.41557e-02 DIIS
   @ROHF iter  52:    -0.81331391614163    1.26210e-12   1.41557e-02 DIIS
   @ROHF iter  53:    -0.81331391614288   -1.25322e-12   1.41557e-02 DIIS
   @ROHF iter  54:    -0.81331391614440   -1.51701e-12   1.41557e-02 DIIS
   @ROHF iter  55:    -0.81331391614309    1.31339e-12   1.41557e-02 DIIS
   @ROHF iter  56:    -0.81331391614240    6.86340e-13   1.41557e-02 DIIS
   @ROHF iter  57:    -0.81331391614167    7.29639e-13   1.41557e-02 DIIS
   @ROHF iter  58:    -0.81331391614118    4.94271e-13   1.41557e-02 DIIS
   @ROHF iter  59:    -0.81331391614172   -5.42233e-13   1.41557e-02 DIIS
   @ROHF iter  60:    -0.81331391614361   -1.89160e-12   1.41557e-02 DIIS
   @ROHF iter  61:    -0.81331391614312    4.90941e-13   1.41557e-02 DIIS
   @ROHF iter  62:    -0.81331391614292    1.93623e-13   1.41557e-02 DIIS
   @ROHF iter  63:    -0.81331391614095    1.97020e-12   1.41557e-02 DIIS
   @ROHF iter  64:    -0.81331391614268   -1.72107e-12   1.41557e-02 DIIS
   @ROHF iter  65:    -0.81331391614231    3.69704e-13   1.41557e-02 DIIS
   @ROHF iter  66:    -0.81331391614196    3.47944e-13   1.41557e-02 DIIS
   @ROHF iter  67:    -0.81331391614231   -3.56604e-13   1.41557e-02 DIIS
   @ROHF iter  68:    -0.81331391614077    1.54432e-12   1.41557e-02 DIIS
   @ROHF iter  69:    -0.81331391614218   -1.40754e-12   1.41557e-02 DIIS
   @ROHF iter  70:    -0.81331391614302   -8.41105e-13   1.41557e-02 DIIS
   @ROHF iter  71:    -0.81331391614172    1.29718e-12   1.41557e-02 DIIS
   @ROHF iter  72:    -0.81331391614132    4.06342e-13   1.41557e-02 DIIS
   @ROHF iter  73:    -0.81331391614108    2.32703e-13   1.41557e-02 DIIS
   @ROHF iter  74:    -0.81331391614245   -1.36358e-12   1.41557e-02 DIIS
   @ROHF iter  75:    -0.81331391614073    1.71396e-12   1.41557e-02 DIIS
   @ROHF iter  76:    -0.81331391614239   -1.65379e-12   1.41557e-02 DIIS
   @ROHF iter  77:    -0.81331391614176    6.27942e-13   1.41557e-02 DIIS
   @ROHF iter  78:    -0.81331391614171    4.72955e-14   1.41557e-02 DIIS
   @ROHF iter  79:    -0.81331391614290   -1.18594e-12   1.41557e-02 DIIS
   @ROHF iter  80:    -0.81331391614243    4.68514e-13   1.41557e-02 DIIS
   @ROHF iter  81:    -0.81331391614313   -7.03881e-13   1.41557e-02 DIIS
   @ROHF iter  82:    -0.81331391614303    1.06581e-13   1.41557e-02 DIIS
   @ROHF iter  83:    -0.81331391614185    1.17795e-12   1.41557e-02 DIIS
   @ROHF iter  84:    -0.81331391614402   -2.16960e-12   1.41557e-02 DIIS
   @ROHF iter  85:    -0.81331391614309    9.22817e-13   1.41557e-02 DIIS
   @ROHF iter  86:    -0.81331391614284    2.53353e-13   1.41557e-02 DIIS
   @ROHF iter  87:    -0.81331391614234    5.00266e-13   1.41557e-02 DIIS
   @ROHF iter  88:    -0.81331391614106    1.28098e-12   1.41557e-02 DIIS
   @ROHF iter  89:    -0.81331391614072    3.38618e-13   1.41557e-02 DIIS
   @ROHF iter  90:    -0.81331391614409   -3.36708e-12   1.41557e-02 DIIS
   @ROHF iter  91:    -0.81331391614210    1.99218e-12   1.41557e-02 DIIS
   @ROHF iter  92:    -0.81331391614264   -5.47562e-13   1.41557e-02 DIIS
   @ROHF iter  93:    -0.81331391614246    1.79856e-13   1.41557e-02 DIIS
   @ROHF iter  94:    -0.81331391614326   -7.96474e-13   1.41557e-02 DIIS
   @ROHF iter  95:    -0.81331391614549   -2.23022e-12   1.41557e-02 DIIS
   @ROHF iter  96:    -0.81331391614341    2.08367e-12   1.41557e-02 DIIS
   @ROHF iter  97:    -0.81331391614133    2.07323e-12   1.41557e-02 DIIS
   @ROHF iter  98:    -0.81331391614135   -1.75415e-14   1.41557e-02 DIIS
   @ROHF iter  99:    -0.81331391614122    1.27010e-13   1.41557e-02 DIIS
   @ROHF iter 100:    -0.81331391614400   -2.77622e-12   1.41557e-02 DIIS

PsiException: Could not converge SCF iterations in 100 iterations.

  Energy and/or wave function did not converge, but proceeding anyway.


  @ROHF Final Energy:    -0.81331391614400

   => Energetics <=

    Nuclear Repulsion Energy =              0.5121069780677421
    One-Electron Energy =                  -2.4334739832061794
    Two-Electron Energy =                   1.1080530889944371
    Total Energy =                         -0.8133139161440002




         ---------------------------------------------------------
                          Configuration Interaction
                            (a 'D E T C I' module)

                 C. David Sherrill, Daniel G. A. Smith, and
                              Matt L. Leininger
         ---------------------------------------------------------

     Iter   Root       Total Energy       Delta E      C RMS

   @CI  0:     0     -1.399986986556   -1.9121E+00   1.2885E-15  
    Warning: Norm of correction (root 0) is < 1.0E-13
   @CI  1:     0     -1.399986986556   0.0000E+00   2.2436E-15 c

   ==> Energetics <==

    SCF energy =           -0.813313916144000
    Total CI energy =      -1.399986986555533

   ==> FCI root 0 information <==

    FCI Root 0 energy =    -1.399986986555533

   The 9 most important determinants:

    *   1    0.557034  (    1,    0)  1AX 3AA 
    *   2    0.524980  (    0,    2)  1AA 2AA 3AB 
    *   3   -0.520979  (    2,    1)  2AX 3AA 
    *   4    0.303121  (    1,    1)  1AA 2AB 3AA 
    *   5   -0.221859  (    2,    0)  1AB 2AA 3AA 
    *   6   -0.031303  (    0,    0)  1AX 2AA 
    *   7   -0.024537  (    1,    2)  1AA 3AX 
    *   8    0.002035  (    0,    1)  1AA 2AX 


We suspect an issue with the Intel compilers used to build the Psi we put on conda.

Turning soscf on should get you around this error.

1 Like

I’d also suggest an input more like the below:

  • energy("fci") will run its own SCF and use detci
  • the e/d/r_convergence will be set a little tighter by default for correlated calcs
molecule mol {
    H 0.0        3.1000000000000014 0.0
    H 0.0 0.0 0.0
    H        2.6846787517317607        1.5500000000000012 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
#    fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf true
    scf_type pk
    maxiter 100
    num_amps_print 1e6
#    r_convergence 1e-6
#    d_convergence 1e-6
#    e_convergence 1e-6
#    ints_tolerance EQUALITY_TOLERANCE
#    damping_percentage 0
}
#set qc_module detci

#hf_energy= energy('scf')
fci_energy = energy('fci')

Thanks for all those helpful comments. Turning on the soscf indeed solved the problem for this particular geometry. However, it seem to be not robust enough to give correct result for some other geometries. Below are two examples:

Output 1

    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.5 release

                         Git: Rev {HEAD} e9f4d6d 

  ==> Input File <==

--------------------------------------------------------------------------

molecule mol {
    H 0.0        3.0000000000000009 0.0
    H 0.0 0.0 0.0
    H        2.5980762113533165        1.5000000000000009 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
   #fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf true
    scf_type pk
    maxiter 100
    num_amps_print 1e6
   # r_convergence 1e-6
   # d_convergence 1e-6
   # e_convergence 1e-6
   # ints_tolerance EQUALITY_TOLERANCE
   # damping_percentage 0
}
# set qc_module detci

fci_energy = energy('fci')


         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H           -0.000000000000     0.000000000000     1.732050807569     1.007825032230
         H           -1.500000000000     0.000000000000    -0.866025403784     1.007825032230
         H            1.500000000000     0.000000000000    -0.866025403784     1.007825032230

  Running in c1 symmetry.


  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter SAD:    -0.81905128680088   -8.19051e-01   0.00000e+00 
   @ROHF iter   1:    -1.03377420625418   -2.14723e-01   3.96859e-02 DIIS
   @ROHF iter   2:    -0.84274319120064    1.91031e-01   2.68562e-02 DIIS
   @ROHF iter   3:    -0.90511173363194   -6.23685e-02   4.79707e-02 DIIS
   @ROHF iter   4:    -0.89551979494072    9.59194e-03   4.62542e-02 DIIS
   @ROHF iter   5:    -1.10124303503219   -2.05723e-01   2.61728e-02 DIIS
   @ROHF iter   6:    -1.07792543029941    2.33176e-02   3.66058e-02 DIIS
   @ROHF iter   7:    -1.01550374807050    6.24217e-02   4.96717e-02 DIIS
   @ROHF iter   8:    -1.02482525806459   -9.32151e-03   4.84216e-02 DIIS
   @ROHF iter   9:    -0.93276769924333    9.20576e-02   5.10174e-02 DIIS
   @ROHF iter  10:    -0.83439668908153    9.83710e-02   3.72907e-02 DIIS
   @ROHF iter  11:    -0.85665689462484   -2.22602e-02   3.66223e-02 DIIS
   @ROHF iter  12:    -0.81399068063143    4.26662e-02   8.77909e-03 DIIS
   @ROHF iter  13:    -0.80904695553229    4.94373e-03   6.67643e-03 SOSCF, nmicro=4
   @ROHF iter  14:    -0.80805650880347    9.90447e-04   2.45908e-03 SOSCF, nmicro=4
   @ROHF iter  15:    -0.80789628176372    1.60227e-04   1.60445e-04 SOSCF, nmicro=4
   @ROHF iter  16:    -0.80789708614447   -8.04381e-07   1.40688e-06 SOSCF, nmicro=4
   @ROHF iter  17:    -0.80789708610573    3.87468e-11   3.68006e-11 SOSCF, nmicro=conv
  Energy and wave function converged.



  @ROHF Final Energy:    -0.80789708610573

   => Energetics <=

    Nuclear Repulsion Energy =              0.5291772106700000
    One-Electron Energy =                  -2.4650068875526028
    Two-Electron Energy =                   1.1279325907768751
    Total Energy =                         -0.8078970861057277

Output 2

    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.5 release

                         Git: Rev {HEAD} e9f4d6d 


  ==> Input File <==

--------------------------------------------------------------------------

molecule mol {
    H 0.0        3.0000000000 0.0
    H 0.0 0.0 0.0
    H        2.5980762113533165        1.500000000000 0.0
symmetry c1
}


set reference rohf

# Set global parameters of calculation.
set globals {
    print_MOs true
    basis sto-3g
    freeze_core false
   #fail_on_maxiter False
    df_scf_guess false
    opdm true
    tpdm true
    soscf true
    scf_type pk
    maxiter 100
    num_amps_print 1e6
   # r_convergence 1e-6
   # d_convergence 1e-6
   # e_convergence 1e-6
   # ints_tolerance EQUALITY_TOLERANCE
   # damping_percentage 0
}
# set qc_module detci

fci_energy = energy('fci')
--------------------------------------------------------------------------



         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H           -0.000000000000     0.000000000000     1.732050807569     1.007825032230
         H           -1.500000000000     0.000000000000    -0.866025403784     1.007825032230
         H            1.500000000000     0.000000000000    -0.866025403784     1.007825032230

  Running in c1 symmetry.


  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter SAD:    -0.81905128680088   -8.19051e-01   0.00000e+00 
   @ROHF iter   1:    -1.03377420625418   -2.14723e-01   4.03667e-02 DIIS
   @ROHF iter   2:    -0.84108463008364    1.92690e-01   2.60508e-02 DIIS
   @ROHF iter   3:    -0.89737454599494   -5.62899e-02   4.68768e-02 DIIS
   @ROHF iter   4:    -0.88507393251099    1.23006e-02   4.44585e-02 DIIS
   @ROHF iter   5:    -1.10832234736944   -2.23248e-01   2.15614e-02 DIIS
   @ROHF iter   6:    -1.11792886346888   -9.60652e-03   1.26774e-02 DIIS
   @ROHF iter   7:    -1.11857784341106   -6.48980e-04   1.18244e-02 DIIS
   @ROHF iter   8:    -1.12194559548220   -3.36775e-03   4.87990e-03 DIIS
   @ROHF iter   9:    -1.12258007362880   -6.34478e-04   1.28959e-03 SOSCF, nmicro=3
   @ROHF iter  10:    -1.12262833029823   -4.82567e-05   7.74474e-06 SOSCF, nmicro=3
   @ROHF iter  11:    -1.12262833206677   -1.76854e-09   7.00835e-09 SOSCF, nmicro=conv
  Energy and wave function converged.

  @ROHF Final Energy:    -1.12262833206677

   => Energetics <=

    Nuclear Repulsion Energy =              0.5291772106700000
    One-Electron Energy =                  -2.4772421024126152
    Two-Electron Energy =                   0.8254365596758491
    Total Energy =                         -1.1226283320667663

Computation Completed

This is a classical case of multiple local minima. HF is really bad at stretched geometries, since the system becomes strongly correlated and thereby affords many local minima for HF. What you probably want to do is run HF with many initial guesses in the hopes of locating the lowest solution, see e.g. doi:10.1063/5.0080817. Another possibility would be to run FCI and use its natural orbitals as the mean-field state; this will give you an unambiguous starting point.

Thanks! This is very helpful. I’m pretty interested in using FCI result as the starting point. But I’m a little confused how to set this up. Can you give me a code snippet or point me to the correct direction?

Here’s an example of using various guesses

molecule {
symmetry c1
H 0.0 3.1000000000000014 0.0
H 0.0 0.0 0.0
H 2.6846787517317607 1.5500000000000012 0.0
}

set reference rohf
set basis sto-3g

# turn off density fitting which is not useful for sto-3g
set scf_type pk
set df_scf_guess false

# turn off crash for demonstration
set fail_on_maxiter false

# core guess
set guess core
energy('scf')

# sad guess (default)
set guess sad
energy('scf')

# huckel guess
set guess huckel
energy('scf')

# gwh guess
set guess gwh
energy('scf')

# sap guess; only available in >= 1.4
set guess sap
energy('scf')

Running this in Psi4 1.4.something (snapshot compiled in March)

  • core: last energy -0.81159653616247, no convergence
  • SAD: converged to energy -1.04864710192611
  • Huckel: converged to energy -1.11610323722485
  • GWH: converged to energy -1.04864710193125
  • SAP: converged to energy -1.04864710193135

As you can see, the Huckel guess leads to the lowest energy.

Repeating the calculation in Psi4 1.3.2, I get

  • core: converged to -1.11610323722483
  • SAD: last energy -1.09227169739243, no convergence
  • Huckel: converged to energy -1.04864710192493
  • GWH: converged to energy -1.11610323722483

and so the behavior is completely different, confirming that this is just a really hard system.

Aha! Thanks for the quick reply. This is very interesting indeed. I understand giving a good guess at low-cost is difficult. Even though some of the guesses above gave the correct result, it is still possible none of the above guesses give correct result.

However, this means it does not quite solve my problem. To conclude this tread, here is what I need:

Suppose currently I have a small system and have the luxury of using FCI. For the benchmark suppose I want to reliably get the ROHF result, how do I set up to use the FCI as initial guess as you suggested? I don’t see a particular keyword for that. Do I need to customize and use READ? Is there a working code that can do this?

This should do what you want

molecule {
symmetry c1
H 0.0 3.1000000000000014 0.0
H 0.0 0.0 0.0
H 2.6846787517317607 1.5500000000000012 0.0
}

# Find FCI natural orbitals
set reference rohf
set basis sto-3g
set guess huckel
set scf_type pk
set df_scf_guess false
set soscf true
set nat_orbs true
energy('fci')

# Read in the natural orbitals as guess to ROHF
set guess read
energy('scf')

but the problem is again that the convergence can be really unpredictable: the calculation started from the FCI natural orbitals looks like it may have a bit sketchy convergence

   @ROHF iter   0:    -0.81047314720426   -8.10473e-01   1.01928e-15 
   @ROHF iter   1:    -1.11610323722485   -3.05630e-01   8.54385e-15 DIIS
   @ROHF iter   2:    -1.11610323722485    8.88178e-16   1.99289e-14 DIIS
  Energy and wave function converged.

It does land on the lowest solution found above.

Hi @susilehtola , thanks for the input. I tried the input and it works for the example shown above, but it does not get me very far, which really confuses me. I use the input format you gave me, with a similar geometry, but it still generates wrong result. Any insight why that might be the case? Does it mean even a FCI guess can’t guarantee the convergence to the ROHF ground state?

I also separately tried two other initial guesses that gave me the following:

Core: last energy -0.80414064804245, on convergence
SAD: converged to energy -1.11029339017942

FCI: -1.399913785979888

Below is the input/output file.

    -----------------------------------------------------------------------
          Psi4: An Open-Source Ab Initio Electronic Structure Package
                               Psi4 1.5 release

                         Git: Rev {HEAD} e9f4d6d 

  ==> Input File <==

--------------------------------------------------------------------------
molecule {
symmetry c1
    H 0.0        3.2000000000000006 0.0
    H 0.0 0.0 0.0
    H        2.7712812921102041        1.6000000000000008 0.0
}

# Find FCI natural orbitals
set reference rohf
set basis sto-3g
set guess huckel
set scf_type pk
set df_scf_guess false
set soscf true
set nat_orbs true
energy('fci')

# Read in the natural orbitals as guess to ROHF
set guess read
energy('scf')--------------------------------------------------------------------------

Scratch directory: /tmp/

*** tstart() called on Kee-Macbook-Pro-16.local
*** at Tue May  3 22:31:55 2022

   => Loading Basis Set <=

    Name: STO-3G
    Role: ORBITAL
    Keyword: BASIS
    atoms 1-3 entry H          line    19 file /Users/qingfengwang/opt/anaconda3/envs/p4env/share/psi4/basis/sto-3g.gbs 


         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H            0.000000000000     0.000000000000     1.847520861407     1.007825032230
         H           -1.600000000000     0.000000000000    -0.923760430703     1.007825032230
         H            1.600000000000     0.000000000000    -0.923760430703     1.007825032230

  Running in c1 symmetry.

  Rotational constants: A =      3.26694  B =      3.26694  C =      1.63347 [cm^-1]
  Rotational constants: A =  97940.45048  B =  97940.45048  C =  48970.22524 [MHz]
  Nuclear repulsion =    0.496103635003125

  Charge       = 0
  Multiplicity = 2
  Electrons    = 3
  Nalpha       = 2
  Nbeta        = 1

  ==> Algorithm <==

  SCF Algorithm Type is PK.
  DIIS enabled.
  MOM disabled.
  Fractional occupation disabled.
  Guess Type is HUCKEL.
  Energy threshold   = 1.00e-08
  Density threshold  = 1.00e-08
  Integral threshold = 1.00e-12

  ==> Primary Basis <==

  Basis Set: STO-3G
    Blend: STO-3G
    Number of shells: 3
    Number of basis functions: 3
    Number of Cartesian functions: 3
    Spherical Harmonics?: true
    Max angular momentum: 0

  ==> Integral Setup <==

  Using in-core PK algorithm.
   Calculation information:
      Number of atoms:                   3
      Number of AO shells:               3
      Number of primitives:              9
      Number of atomic orbitals:         3
      Number of basis functions:         3

      Integral cutoff                 1.00e-12
      Number of threads:                 1

  Performing in-core PK
  Using 42 doubles for integral storage.
  We computed 21 shell quartets total.
  Whereas there are 21 unique shell quartets.

  ==> DiskJK: Disk-Based J/K Matrices <==

    J tasked:                  Yes
    K tasked:                  Yes
    wK tasked:                  No
    Memory [MiB]:              375
    Schwarz Cutoff:          1E-12

    OpenMP threads:              1

  Minimum eigenvalue in the overlap matrix is 9.8813656671E-01.
  Reciprocal condition number of the overlap matrix is 9.6523457474E-01.
    Using symmetric orthogonalization.

  ==> Pre-Iterations <==

  SCF Guess: Huckel guess via on-the-fly atomic UHF (doi:10.1021/acs.jctc.8b01089).

   -------------------------------------------------------
    Irrep   Nso     Nmo     Nalpha   Nbeta   Ndocc  Nsocc
   -------------------------------------------------------
     A          3       3       2       1       1       1
   -------------------------------------------------------
    Total       3       3       2       1       1       1
   -------------------------------------------------------

  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter   1:    -1.01645810957191   -1.01646e+00   4.72321e-02 DIIS
   @ROHF iter   2:    -0.80896680156876    2.07491e-01   1.56561e-02 DIIS
   @ROHF iter   3:    -0.80778037437229    1.18643e-03   2.04625e-02 DIIS
   @ROHF iter   4:    -0.80524884308400    2.53153e-03   1.83448e-02 DIIS
   @ROHF iter   5:    -0.80661608907971   -1.36725e-03   1.37592e-02 DIIS
   @ROHF iter   6:    -0.81275015803598   -6.13407e-03   1.88275e-02 DIIS
   @ROHF iter   7:    -1.07908774290014   -2.66338e-01   2.92358e-02 DIIS
   @ROHF iter   8:    -0.80882080895932    2.70267e-01   1.55556e-02 DIIS
   @ROHF iter   9:    -0.80394478900721    4.87602e-03   1.11884e-02 DIIS
   @ROHF iter  10:    -0.80765576113253   -3.71097e-03   1.45895e-02 DIIS
   @ROHF iter  11:    -0.80667249412817    9.83267e-04   1.37736e-02 DIIS
   @ROHF iter  12:    -0.79793284998879    8.73964e-03   9.56295e-03 DIIS
   @ROHF iter  13:    -0.79959542395487   -1.66257e-03   1.21803e-02 SOSCF, nmicro=3
   @ROHF iter  14:    -0.79530002274239    4.29540e-03   2.48841e-04 DIIS
   @ROHF iter  15:    -0.80197110628307   -6.67108e-03   8.69827e-03 SOSCF, nmicro=2
   @ROHF iter  16:    -0.79891425384340    3.05685e-03   1.05127e-04 SOSCF, nmicro=4
   @ROHF iter  17:    -0.79891469696071   -4.43117e-07   1.83838e-06 SOSCF, nmicro=4
   @ROHF iter  18:    -0.79891469688663    7.40814e-11   1.00127e-10 SOSCF, nmicro=conv
  Energy and wave function converged.


  ==> Post-Iterations <==

    Orbital Energies [Eh]
    ---------------------

    Doubly Occupied:                                                      

       1A     -0.300720  

    Singly Occupied:                                                      

       2A     -0.018695  

    Virtual:                                                              

       3A      0.080952  

    Final Occupation by Irrep:
              A 
    DOCC [     1 ]
    SOCC [     1 ]

  @ROHF Final Energy:    -0.79891469688663

   => Energetics <=

    Nuclear Repulsion Energy =              0.4961036350031248
    One-Electron Energy =                  -2.3970115343338381
    Two-Electron Energy =                   1.1019932024440817
    Total Energy =                         -0.7989146968866314

Computation Completed


Properties will be evaluated at   0.000000,   0.000000,   0.000000 [a0]

Properties computed using the SCF density matrix

  Nuclear Dipole Moment: [e a0]
     X:    -0.0000      Y:     0.0000      Z:     0.0000

  Electronic Dipole Moment: [e a0]
     X:     4.5334      Y:     0.0000      Z:     2.6173

  Dipole Moment: [e a0]
     X:     4.5334      Y:     0.0000      Z:     2.6173     Total:     5.2347

  Dipole Moment: [D]
     X:    11.5226      Y:     0.0000      Z:     6.6526     Total:    13.3052


*** tstop() called on Kee-Macbook-Pro-16.local at Tue May  3 22:31:56 2022
Module time:
	user time   =       0.16 seconds =       0.00 minutes
	system time =       0.01 seconds =       0.00 minutes
	total time  =          1 seconds =       0.02 minutes
Total time:
	user time   =       0.16 seconds =       0.00 minutes
	system time =       0.01 seconds =       0.00 minutes
	total time  =          1 seconds =       0.02 minutes
 MINTS: Wrapper to libmints.
   by Justin Turney

   Calculation information:
      Number of threads:                 1
      Number of atoms:                   3
      Number of AO shells:               3
      Number of SO shells:               3
      Number of primitives:              9
      Number of atomic orbitals:         3
      Number of basis functions:         3

      Number of irreps:                  1
      Integral cutoff                 1.00e-12
      Number of functions per irrep: [   3 ]

 OEINTS: Overlap, kinetic, potential, dipole, and quadrupole integrals
         stored in file 35.

      Computing two-electron integrals...done
      Computed 21 non-zero two-electron integrals.
        Stored in file 33.


         ---------------------------------------------------------
                          Configuration Interaction
                            (a 'D E T C I' module)

                 C. David Sherrill, Daniel G. A. Smith, and
                              Matt L. Leininger
         ---------------------------------------------------------

   ==> Starting CI iterations <==

    H0 Block Eigenvalue =  -1.39991379

    Simultaneous Expansion Method (Block Davidson Method)
    Using 1 initial trial vectors

     Iter   Root       Total Energy       Delta E      C RMS

   @CI  0:     0     -1.399913785980   -1.8960E+00   8.6785E-16  
    Warning: Norm of correction (root 0) is < 1.0E-13
   @CI  1:     0     -1.399913785980   0.0000E+00   1.0010E-15 c

   Computing CI Natural Orbitals

   ==> Energetics <==

    SCF energy =           -0.798914696886631
    Total CI energy =      -1.399913785979888

   ==> FCI root 0 information <==

    FCI Root 0 energy =    -1.399913785979888

   Active Space Natural occupation numbers:

         A   1.052998         A   1.000000         A   0.947002

   The 9 most important determinants:

    *   1   -0.815502  (    2,    0)  1AB 2AA 3AA 
    *   2   -0.407800  (    1,    1)  1AA 2AB 3AA 
    *   3    0.407702  (    0,    2)  1AA 2AA 3AB 
    *   4   -0.049020  (    1,    0)  1AX 3AA 
    *   5    0.005740  (    2,    1)  2AX 3AA 
    *   6   -0.000000  (    0,    1)  1AA 2AX 
    *   7    0.000000  (    1,    2)  1AA 3AX 
    *   8   -0.000000  (    0,    0)  1AX 2AA 

		 "A good bug is a dead bug" 

			 - Starship Troopers

		 "I didn't write FORTRAN.  That's the problem."

			 - Edward Valeev

Scratch directory: /tmp/

*** tstart() called on Kee-Macbook-Pro-16.local
*** at Tue May  3 22:31:56 2022

   => Loading Basis Set <=

    Name: STO-3G
    Role: ORBITAL
    Keyword: BASIS
    atoms 1-3 entry H          line    19 file /Users/qingfengwang/opt/anaconda3/envs/p4env/share/psi4/basis/sto-3g.gbs 


         ---------------------------------------------------------
                                   SCF
               by Justin Turney, Rob Parrish, Andy Simmonett
                          and Daniel G. A. Smith
                             ROHF Reference
                        1 Threads,    500 MiB Core
         ---------------------------------------------------------

  ==> Geometry <==

    Molecular point group: c1
    Full point group: D3h

    Geometry (in Angstrom), charge = 0, multiplicity = 2:

       Center              X                  Y                   Z               Mass       
    ------------   -----------------  -----------------  -----------------  -----------------
         H            0.000000000000     0.000000000000     1.847520861407     1.007825032230
         H           -1.600000000000     0.000000000000    -0.923760430703     1.007825032230
         H            1.600000000000     0.000000000000    -0.923760430703     1.007825032230

  Running in c1 symmetry.

  Rotational constants: A =      3.26694  B =      3.26694  C =      1.63347 [cm^-1]
  Rotational constants: A =  97940.45048  B =  97940.45048  C =  48970.22524 [MHz]
  Nuclear repulsion =    0.496103635003125

  Charge       = 0
  Multiplicity = 2
  Electrons    = 3
  Nalpha       = 2
  Nbeta        = 1

  ==> Algorithm <==

  SCF Algorithm Type is PK.
  DIIS enabled.
  MOM disabled.
  Fractional occupation disabled.
  Guess Type is READ.
  Energy threshold   = 1.00e-06
  Density threshold  = 1.00e-06
  Integral threshold = 1.00e-12

  ==> Primary Basis <==

  Basis Set: STO-3G
    Blend: STO-3G
    Number of shells: 3
    Number of basis functions: 3
    Number of Cartesian functions: 3
    Spherical Harmonics?: true
    Max angular momentum: 0

   => Loading Basis Set <=

    Name: STO-3G
    Role: ORBITAL
    Keyword: BASIS
    atoms 1-3 entry H          line    19 file /Users/qingfengwang/opt/anaconda3/envs/p4env/share/psi4/basis/sto-3g.gbs 

  Reading orbitals from file /tmp/1_r=3.2.input.default.11497.180.npy, no projection.

  ==> Integral Setup <==

  Using in-core PK algorithm.
   Calculation information:
      Number of atoms:                   3
      Number of AO shells:               3
      Number of primitives:              9
      Number of atomic orbitals:         3
      Number of basis functions:         3

      Integral cutoff                 1.00e-12
      Number of threads:                 1

  Performing in-core PK
  Using 42 doubles for integral storage.
  We computed 21 shell quartets total.
  Whereas there are 21 unique shell quartets.

  ==> DiskJK: Disk-Based J/K Matrices <==

    J tasked:                  Yes
    K tasked:                  Yes
    wK tasked:                  No
    Memory [MiB]:              375
    Schwarz Cutoff:          1E-12

    OpenMP threads:              1

  Minimum eigenvalue in the overlap matrix is 9.8813656671E-01.
  Reciprocal condition number of the overlap matrix is 9.6523457474E-01.
    Using symmetric orthogonalization.

  ==> Pre-Iterations <==

  SCF Guess: Orbitals guess was supplied from a previous computation.

   -------------------------------------------------------
    Irrep   Nso     Nmo     Nalpha   Nbeta   Ndocc  Nsocc
   -------------------------------------------------------
     A          3       3       2       1       1       1
   -------------------------------------------------------
    Total       3       3       2       1       1       1
   -------------------------------------------------------

  ==> Iterations <==

                        Total Energy        Delta E     RMS |[F,P]|

   @ROHF iter   0:    -0.79646150952498   -7.96462e-01   6.35853e-03 
   @ROHF iter   1:    -0.80115077426506   -4.68926e-03   7.44594e-03 DIIS
   @ROHF iter   2:    -0.79784088930282    3.30988e-03   9.38835e-03 DIIS
   @ROHF iter   3:    -0.84719987358048   -4.93590e-02   3.89501e-02 SOSCF, nmicro=2
   @ROHF iter   4:    -0.80103402710586    4.61658e-02   1.40344e-02 DIIS
   @ROHF iter   5:    -0.80223181990670   -1.19779e-03   9.05636e-03 DIIS
   @ROHF iter   6:    -0.79585436204199    6.37746e-03   4.39281e-03 SOSCF, nmicro=2
   @ROHF iter   7:    -0.79530204454024    5.52318e-04   1.03171e-05 SOSCF, nmicro=2
   @ROHF iter   8:    -0.79530204149895    3.04129e-09   2.75288e-10 SOSCF, nmicro=conv
  Energy and wave function converged.


  ==> Post-Iterations <==

    Orbital Energies [Eh]
    ---------------------

    Doubly Occupied:                                                      

       1A     -0.238741  

    Singly Occupied:                                                      

       2A     -0.142792  

    Virtual:                                                              

       3A      0.143042  

    Final Occupation by Irrep:
              A 
    DOCC [     1 ]
    SOCC [     1 ]

  @ROHF Final Energy:    -0.79530204149895

   => Energetics <=

    Nuclear Repulsion Energy =              0.4961036350031248
    One-Electron Energy =                  -2.3971528878894244
    Two-Electron Energy =                   1.1057472113873499
    Total Energy =                         -0.7953020414989496

Computation Completed


Properties will be evaluated at   0.000000,   0.000000,   0.000000 [a0]

Properties computed using the SCF density matrix

  Nuclear Dipole Moment: [e a0]
     X:     0.0000      Y:     0.0000      Z:     0.0000

  Electronic Dipole Moment: [e a0]
     X:    -4.5228      Y:     0.0000      Z:    -2.6112

  Dipole Moment: [e a0]
     X:    -4.5228      Y:     0.0000      Z:    -2.6112     Total:     5.2224

  Dipole Moment: [D]
     X:   -11.4957      Y:     0.0000      Z:    -6.6371     Total:    13.2741


*** tstop() called on Kee-Macbook-Pro-16.local at Tue May  3 22:31:56 2022
Module time:
	user time   =       0.08 seconds =       0.00 minutes
	system time =       0.00 seconds =       0.00 minutes
	total time  =          0 seconds =       0.00 minutes
Total time:
	user time   =       0.25 seconds =       0.00 minutes
	system time =       0.02 seconds =       0.00 minutes
	total time  =          1 seconds =       0.02 minutes

    Psi4 stopped on: Tuesday, 03 May 2022 10:31PM
    Psi4 wall time for execution: 0:00:00.41

*** Psi4 exiting successfully. Buy a developer a beer!



Well, H3 in minimal basis is a 3-orbital problem, so it should be even feasible to do a brute force search of the orbital space. One just needs to sample the rotation angles to go through the whole solution space. I think that is the ultimate solution for this problem.

1 Like

Thanks for the suggestion.
For my purpose, since I want the ROHF energy in a grid of geometry, I can take advantage of that and start with a geometry that gives correct ROHF energy, and feed the C matrix to the next neighbor geometry. This is sort of like adiabatic approach. It works for me and gives me correct and smooth ROHF energy alone the dissociation curve.
I haven’t tested the brute force method you mentioned, but I expect it will work as well.