Move to LibECPInt (#2135)

* Add beginnings of LibECPInt external build infrastructure

* Add (broken) libECP energy calls

* Apply fix to code suggested by RAS, more debug info

* Fix order order of l and n

* Apply suggestions from @robashaw, which fix the energies

* Remove old ECP implementation

* Implement analytic ECP gradients

* Make ECP gradients a little more efficient

* Initial ECP skeleton second derivative implementation

* Correct bookkeeping errors, and second derivative call

* Rebase to L2 codebase

* SCF ECP test and auto libecpint build working

* Reimplement spherical transform

* Implement RHF and UHF ECP Hessians

* ecpint testing marks and optional

* ecp allow disabled, add ecosys gha

* ci 1

* ci 2

* ci 3

* relax cppe+ecp to rest of file. remove extra ecp build

Co-authored-by: Lori A. Burns <[email protected]>

.github/workflows/ecosystem.yml

@@ -118,6 +118,7 @@ jobs:
  # qc opt'l
  - psi4/label/dev::ambit
  - psi4/label/dev::chemps2
+ - psi4/label/dev::libecpint
  - psi4/label/dev::dkh
  - psi4/label/dev::gdma
  - psi4/label/dev::pcmsolver
@@ -147,6 +148,7 @@ jobs:
  # nyi
  sed -i "s;- psi4/label/dev::ambit;;g" build.yaml
  sed -i "s;- psi4/label/dev::chemps2;;g" build.yaml
+ sed -i "s;- psi4/label/dev::libecpint;;g" build.yaml
  sed -i "s;- psi4/label/dev::dkh;;g" build.yaml
  sed -i "s;- psi4/label/dev::gdma;;g" build.yaml
  sed -i "s;- psi4/label/dev::pcmsolver;;g" build.yaml
@@ -180,6 +182,7 @@ jobs:
  - psi4/label/dev::ambit
  - psi4/label/dev::chemps2
  - psi4/label/dev::dkh
+ - psi4/label/dev::libecpint
  - psi4/label/dev::gdma
  - psi4/label/dev::pcmsolver
  - psi4/label/dev::simint
@@ -249,6 +252,7 @@ jobs:
  - libxc
  # qc req'd from buildtime
  #- psi4/label/dev::ambit
+ #- conda-forge::libecpint
  #- psi4/label/dev::chemps2
  #- psi4/label/dev::dkh
  #- psi4/label/dev::gdma
@@ -328,6 +332,8 @@ jobs:
  -D CMAKE_INSIST_FIND_PACKAGE_CheMPS2=ON \
  -D ENABLE_dkh=ON \
  -D CMAKE_INSIST_FIND_PACKAGE_dkh=ON \
+ -D ENABLE_ecpint=ON \
+ -D CMAKE_INSIST_FIND_PACKAGE_ecpint=ON \
  -D ENABLE_gdma=ON \
  -D CMAKE_INSIST_FIND_PACKAGE_gdma=ON \
  -D ENABLE_PCMSolver=ON \
@@ -366,6 +372,8 @@ jobs:
  -D CMAKE_INSIST_FIND_PACKAGE_CheMPS2=ON \
  -D ENABLE_dkh=ON \
  -D CMAKE_INSIST_FIND_PACKAGE_dkh=ON \
+ -D ENABLE_ecpint=ON \
+ -D CMAKE_INSIST_FIND_PACKAGE_ecpint=ON \
  -D ENABLE_gdma=ON \
  -D CMAKE_INSIST_FIND_PACKAGE_gdma=ON \
  -D ENABLE_PCMSolver=ON \

CMakeLists.txt

@@ -100,7 +100,7 @@ option(ENABLE_CheMPS2 "Enables CheMPS2 for DMRG (requires HDF5)" OFF)
 option(ENABLE_cppe "Enables CPPE for Polarizable Embedding" OFF)
 option(ENABLE_adcc "Enables adcc for algebraic-diagrammatic construction methods (can also be added at runtime)" OFF)
 option(ENABLE_dkh "Enables DKH integrals (requires Fortran)" OFF)
-option(ENABLE_ecpint "Enables Shaw's libecpint project needed for effective core potentials" OFF)
+option(ENABLE_ecpint "Enables libecpint for effective core potentials (ECP)" OFF)
 option(ENABLE_libefp "Enables LIBEFP and PylibEFP for fragments" OFF)
 option(ENABLE_erd "Enables use of ERD instead of Libint (requires Fortran)" OFF)
 option(ENABLE_simint "Enables use of SIMINT two-electron integral library" OFF)

doc/sphinxman/source/build_planning.rst

@@ -280,7 +280,9 @@ are available pre-built from conda.
 
 .. * :ref:`erd <cmake:erd>` |w---w| :ref:`[what is this?] <sec:erd>` `[min version] <https://github.com/psi4/psi4/blob/master/external/upstream/erd/CMakeLists.txt#L2>`_
 
- * :ref:`Fortran Compiler <cmake:fortran>`
+.. * :ref:`Fortran Compiler <cmake:fortran>`
+
+* :ref:`ecpint <cmake:ecpint>` |w---w| :ref:`[what is ecpint?] <sec:ecpint>` :source:`[ecpint min version] <external/upstream/ecpint/CMakeLists.txt/#L2>`
 
 * :ref:`dkh <cmake:dkh>` |w---w| :ref:`[what is dkh?] <sec:dkh>` :source:`[dkh min version] <external/upstream/dkh/CMakeLists.txt#L2>`
 

external/upstream/ecpint/CMakeLists.txt

@@ -1,5 +1,5 @@
 if(${ENABLE_ecpint})
- find_package(ecpint 1.0.6 CONFIG QUIET)
+ find_package(ecpint 1.0.7 CONFIG QUIET)
 
  if(${ecpint_FOUND})
  get_property(_loc TARGET ECPINT::ecpint PROPERTY LOCATION)
@@ -13,8 +13,8 @@ if(${ENABLE_ecpint})
  include(ExternalProject)
  message(STATUS "Suitable ecpint could not be located, ${Magenta}Building ecpint${ColourReset} instead.")
  ExternalProject_Add(ecpint_external
- URL https://github.com/loriab/libecpint/archive/psiecp.tar.gz
- #UPDATE_COMMAND ""
+ URL https://github.com/robashaw/libecpint/archive/v1.0.7.tar.gz
+ UPDATE_COMMAND ""
  CMAKE_ARGS -DCMAKE_INSTALL_PREFIX=${STAGED_INSTALL_PREFIX}
  -DCMAKE_BUILD_TYPE=${CMAKE_BUILD_TYPE}
  -DCMAKE_C_COMPILER=${CMAKE_C_COMPILER}
@@ -33,7 +33,7 @@ if(${ENABLE_ecpint})
  CMAKE_CACHE_ARGS -DCMAKE_C_FLAGS:STRING=${CMAKE_C_FLAGS}
  -DCMAKE_CXX_FLAGS:STRING=${CMAKE_CXX_FLAGS})
 
- set(ecpint_DIR ${STAGED_INSTALL_PREFIX}/lib/cmake/ecpint CACHE PATH "path to internally built ecpintConfig.cmake" FORCE)
+ set(ecpint_DIR ${STAGED_INSTALL_PREFIX}/lib/cmake/ecpint CACHE PATH "path to internally built ecpint-config.cmake" FORCE)
  endif()
 else()
  add_library(ecpint_external INTERFACE) # dummy

psi4/CMakeLists.txt

@@ -101,7 +101,7 @@ else()
 endif()
 
 if(${ENABLE_ecpint})
- find_package(ecpint 1.0.6 CONFIG REQUIRED)
+ find_package(ecpint 1.0.7 CONFIG REQUIRED)
  get_property(_loc TARGET ECPINT::ecpint PROPERTY LOCATION)
  list(APPEND _addons ${_loc})
  message(STATUS "${Cyan}Using ecpint${ColourReset}: ${_loc} (version ${ecpint_VERSION})")

psi4/driver/procrouting/proc.py

@@ -1451,7 +1451,7 @@ def scf_helper(name, post_scf=True, **kwargs):
  raise ValidationError("""Detected options to both cast and perform a broken symmetry computation""")
 
  if (core.get_option('SCF', 'STABILITY_ANALYSIS') == 'FOLLOW') and (core.get_option('SCF', 'REFERENCE') != 'UHF'):
- raise ValidationError("""Stability analysis root following is only available for UHF""")
+ raise ValidationError(f"""Stability analysis root following is only available for unrestricted Hartree--Fock, not present {core.get_option('SCF', 'REFERENCE')}""")
 
  # broken set-up
  if do_broken:
@@ -2506,45 +2506,6 @@ def run_scf_gradient(name, **kwargs):
 
  grad = core.scfgrad(ref_wfn)
 
- if ref_wfn.basisset().has_ECP():
- core.print_out("\n\n ==> Adding ECP gradient terms (computed numerically) <==\n")
- # Build a map of atom->ECP number
- old_print = ref_wfn.get_print()
- ref_wfn.set_print(0)
- delta = 0.0001
- natom = ref_wfn.molecule().natom()
- mints = core.MintsHelper(ref_wfn)
- ecpgradmat = core.Matrix("ECP Gradient", natom, 3)
- ecpgradmat.zero()
- ecpgrad = np.asarray(ecpgradmat)
- Dmat = ref_wfn.Da_subset("AO")
- Dmat.add(ref_wfn.Db_subset("AO"))
- def displaced_energy(atom, displacement):
- mints.basisset().move_atom(atom, displacement)
- E = Dmat.vector_dot(mints.ao_ecp())
- mints.basisset().move_atom(atom, -1*displacement)
- return E
-
- for atom in range(natom):
- for xyz in range(3):
- transvec = core.Vector3(0.0)
- transvec[xyz] += delta
- # +1 displacement
- Ep1 = displaced_energy(atom, 1*transvec)
- # -1 displacement
- Em1 = displaced_energy(atom, -1*transvec)
- # +2 displacement
- Ep2 = displaced_energy(atom, 2*transvec)
- # -2 displacement
- Em2 = displaced_energy(atom, -2*transvec)
- # Evaluate
- ecpgrad[atom, xyz] = (Em2 + 8*Ep1 - 8*Em1 - Ep2) / (12*delta)
- ecpgradmat.symmetrize_gradient(ref_wfn.molecule())
- ecpgradmat.print_atom_vector()
- grad.add(ecpgradmat)
- grad.print_atom_vector()
- ref_wfn.set_print(old_print)
-
  ref_wfn.set_gradient(grad)
 
  ref_wfn.set_variable("SCF TOTAL GRADIENT", grad) # P::e SCF

psi4/src/CMakeLists.txt

@@ -94,6 +94,17 @@ target_link_libraries(core
  tgt::lapack # OpenMP flags to core.cc
  )
 
+if(TARGET ECPINT::ecpint)
+ target_compile_definitions(core
+ PUBLIC
+ USING_ecpint
+ )
+ target_compile_definitions(l2export
+ PUBLIC
+ USING_ecpint
+ )
+endif()
+
 if(TARGET BrianQC::static_wrapper)
  target_compile_definitions(core
  PRIVATE

psi4/src/export_mints.cc

@@ -141,6 +141,9 @@ std::shared_ptr<BasisSet> construct_basisset_from_pydict(const std::shared_ptr<M
  // basisname is uniform; fill map with key/value (gbs entry) pairs of elements from pybs['shell_map']
  int totalncore = 0;
  if (pybs.contains("ecp_shell_map")) {
+#ifndef USING_ecpint
+ throw PSIEXCEPTION("BasisSet contains ECP shells but libecpint addon not enabled. Re-compile with `-D ENABLE_ecpint=ON`.");
+#endif
  py::list ecpbasisinfo = pybs["ecp_shell_map"].cast<py::list>();
  for (int atom = 0; atom < py::len(ecpbasisinfo); ++atom) {
  std::vector<ShellInfo> vec_shellinfo;
@@ -1430,9 +1433,11 @@ void export_mints(py::module& m) {
  .def("ao_potential", oneelectron_mixed_basis(&MintsHelper::ao_potential), "AO mixed basis potential integrals")
  .def("so_potential", &MintsHelper::so_potential, "SO basis potential integrals",
  "include_perturbations"_a = true)
+#ifdef USING_ecpint
  .def("ao_ecp", oneelectron(&MintsHelper::ao_ecp), "AO basis effective core potential integrals.")
  .def("ao_ecp", oneelectron_mixed_basis(&MintsHelper::ao_ecp), "AO basis effective core potential integrals.")
  .def("so_ecp", &MintsHelper::so_ecp, "SO basis effective core potential integrals.")
+#endif
 
  // One-electron properties and
  .def("ao_pvp", &MintsHelper::ao_pvp, "AO pvp integrals")

psi4/src/psi4/libmints/CMakeLists.txt

@@ -44,9 +44,6 @@ list(APPEND sources
  cdsalclist.cc
  erd_eri.cc
  angularmomentum.cc
- bessel.cc
- gaussquad.cc
- ecpint.cc
  orthog.cc
  )
 
@@ -93,7 +90,6 @@ target_link_libraries(mints
  PUBLIC
  Libint2::cxx
  )
-
 # Conditionally linked-to external projects
 if(TARGET Libint::libint)
  target_link_libraries(mints
@@ -102,6 +98,29 @@ if(TARGET Libint::libint)
  )
 endif()
 
+if(TARGET ECPINT::ecpint)
+ target_sources(mints
+ PRIVATE
+ ecpint.cc
+ )
+ target_compile_definitions(mints
+ PUBLIC
+ USING_ecpint
+ )
+ target_compile_definitions(l2intf
+ PUBLIC
+ USING_ecpint
+ )
+ target_link_libraries(mints
+ PUBLIC
+ ECPINT::ecpint
+ )
+ target_link_libraries(l2intf
+ PUBLIC
+ ECPINT::ecpint
+ )
+endif()
+
 if(TARGET simint::simint)
  target_sources(l2intf
  PRIVATE
@@ -125,13 +144,6 @@ if(TARGET dkh::dkh)
  )
 endif()
 
-if(TARGET ECPINT::ecpint)
- target_link_libraries(mints
- PUBLIC
- ECPINT::ecpint
- )
-endif()
-
 if(TARGET PCMSolver::pcm)
  target_link_libraries(mints
  PUBLIC