PSI4
PSI4 is a C++/Python core that easily interfaces with and is extended by standalone community projects. The major advantage of using PSI4 in tandem with PyOpenCAP is that calculations can be performed in one-shot within the same python script. Since PSI4 allows direct control over data structures such as density matrices, the interface between PSI4 and PyOpenCAP is seamless. Our interface has been tested for the Psi4 dev build, which is available via conda:
conda install -c psi4/label/dev psi4
Step 1: Defining the System object
Molden(recommended)
The best way to construct the System object is to import the geometry
and basis set from a molden file generated by a PSI4. This ensures proper ordering of the
AO basis set.
psi4.molden(wfn, 'molden_in.molden')
molden_dict = {"basis_file":"molden_in.molden","molecule": "molden"}
s = pyopencap.System(molden_dict)
Inline
The molecule and basis set can also be specified inline. The “molecule” keyword must be set to “read”, and then an additional keyword “geometry” must be specified, with a string that contains the geometry in xyz format. The “basis_file” keyword must be set to a path to a basis set file formatted in Psi4 style, which can be downloaded from the MolSSI BSE. Other optional keyword for this section include “bohr_coordinates” and “cart_bf”. Please see the keywords section for more details. It is recommended to check the overlap matrix to ensure that the ordering and normalization matches. Up to G-type functions are supported.
E, wfn = psi4.energy('scf', return_wfn=True)
mints = psi4.core.MintsHelper(wfn.basisset())
S_mat = np.asarray(mints.ao_overlap())
sys_dict = {"geometry": '''N 0 0 1.039
N 0 0 -1.039
X 0 0 0.0''',
"molecule" : "read",
"basis_file":"path/to/basis.bas",
"cart_bf":"d",
"bohr_coordinates:": "true"}
s.check_overlap_mat(S_mat,"psi4")
Step 1: Defining the CAP object
The CAP matrix is computed by the CAP object. The constructor
requires a System object, a dictionary containing the CAP parameters,
and the number of states. An example is provided below. Please see the keywords section for more information on
the CAP parameters.
cap_dict = {"cap_type": "box",
"cap_x":"2.76",
"cap_y":"2.76",
"cap_z":"4.88",
"Radial_precision": "14",
"angular_points": "110"}
nstates = 10
pc = pyopencap.CAP(my_system,cap_dict,nstates)
Step 2: Passing the density matrices
The simplest interface is with the full CI module. One can request one particle densities to be calculated by using the opdm and tdm options:
psi4.set_options({"opdm":True,"num_roots":nstates,"tdm":True,"dipmom":True})
ci_energy, ci_wfn = psi4.energy('FCI', return_wfn=True)
Densities are now available through the get_opdm function. One must be careful to ensure
that the densities are represented in AO basis before passing to PyOpenCAP using the
add_tdm() function:
for i in range(0,nstates):
for j in range(i,nstates):
opdm_mo = ci_wfn.get_opdm(i, j, "SUM", True)
opdm_so = psi4.core.triplet(ci_wfn.Ca(), opdm_mo, ci_wfn.Ca(), False, False, True)
opdm_ao = psi4.core.Matrix(n_bas,n_bas)
opdm_ao.remove_symmetry(opdm_so,so2ao)
pc.add_tdm(opdm_ao.to_array(),i,j,"psi4")
if not i==j:
pc.add_tdm(opdm_ao.to_array().conj().T,j,i,"psi4")
Please see the PSI4 documentation for more details, or our repository for an example.
Note:
The interface with Psi4 is not restricted to FCI. The add_tdm()
function is completely general; it requires only that the densities are in AO basis, and that
the basis set ordering matches the system. An example for ADC is provided in the repository.
Step 3: Computing the CAP matrix
Once all of the densities are loaded, the CAP matrix is computed
using the compute_projected_cap() function. The matrix can be retrieved using the
get_projected_cap() function.
pc.compute_projected_cap()
W_mat=pc.get_projected_cap()
Step 4: Generate and analyze eigenvalue trajectories
H0 and W can be used to construct a CAPHamiltonian object.
In many cases, it can be advantageous to use the
export() function, which generates an OpenCAP
formatted output file, which can be used for later analysis.
from pyopencap.analysis import CAPHamiltonian
CAPH = CAPHamiltonian(H0=H0,W=W_mat)
CAPH.export("output.out")
See the analysis section for more details.