Columbus

COLUMBUS is a collection of programs designed primarily for multi-reference (MR) calculations on electronic ground and excited states of atoms and molecules. Here, we outline the steps of performing CAP-MRCI calculations with Columbus and OpenCAP, though the steps are broadly applicable to any of the methods. These steps have only been tested for the serial version of Columbus.

Step 1: Running MRCI calculation

When running the MRCI calculation, ensure that transition moments between each pair of relevant states are requested. Once the MRCI calculation is finished, navigate to the WORK directory. Assuming one has set up the input properly, the following files will be needed

cid1trfl: files for each pair of states, including state density matrices (i.e. cid1trfl.FROMdrt1.state1TOdrt1.state1)
ciudgsm: file which contains final MRCI energies and convergence information
tranls: file which contains information on active space/frozen orbitals, which is necessary to fully reconstruct density matrices in AO basis

Also, the MO coefficients will be needed, which are located in the MOLDEN directory:

molden_mo_mc.sp: optimized MO coefficients from MCSCF calculation in .molden format

Step 2: Generating human readable density matrix files

The next step is to convert the cid1trfl files into a human readable format. The Columbus utility iwfmt.x can be used for this purpose. We provide a bash script in the main repository utilities/write_iwfmt.bash which can be executed in the WORK directory to generate these files, assuming that the $COLUMBUS environment variable is properly set.

Step 2: Importing the data to PyOpenCAP

System object

To run a PyOpenCAP calculation, the geometry and basis set must be imported into a System object. The constructor takes in a Python dictionary as an argument.

Molden (recommended)

Molden files generated by Columbus contain the geometry and basis set information. To read in from molden, “molden” must be set as the value to the key “molecule”, and the path to the file must be set as the value to the key “basis_file”. Here is an example:

sys_dict = {"molecule": "molden","basis_file": "path/to/file.molden"}
my_system = pyopencap.System(sys_dict)

Inline(not recommended)

The molecule and basis set can also be specified manually. 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. Up to G-type functions are supported.

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"}
my_system = pyopencap.System(sys_dict)

One particle densities/zeroth order Hamiltonian

The CAP matrix is computed by the CAP object. The constructor requires a System, a dictionary containing the CAP parameters, the number of states, and finally the string “openmolcas”, which denotes the ordering of the atomic orbital basis set. 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"}
nstates = 10
pc = pyopencap.CAP(my_system,cap_dict,nstates)

Before we can compute the CAP matrix in the state basis, we must load in the density matrices. Due to the large number of files generated by Columbus, we have provided a colparser utility to manage the data.

A colparser object is instantiated using the tranls file and the MO coefficients:

parser = colparser('data_files/molden_mo_mc.sp', 'data_files/tranls')

The zeroth order Hamiltonian, which is diagonal for MR-CI, can be read in from the ciudgsm file as follows:

H0 = parser.get_H0(filename='data_files/ciudgsm')

Densities are loaded in one at a time using pyopencap.analysis.colparser.sdm_ao() / pyopencap.analysis.colparser.tdm_ao() and add_tdm(). To specify which tdm/sdm to parse, one can use state and optionally DRT indices:

for i in range(0,nstates):
        for j in range(i,nstates):
                if i==j:
                        # Indices start from 0 in pyopencap, but from 1 in Columbus file names
                        dm1_ao = parser.sdm_ao(i+1,data_dir='data_files',DRTn=1)
                        pc.add_tdm(dm1_ao,i,j,'molden')
                else:
                        # Indices start from 0 in pyopencap, but from 1 in Columbus file names
                        dm1_ao = parser.tdm_ao(i+1, j+1,drtFrom=1,drtTo=1,data_dir='data_files')
                        pc.add_tdm(dm1_ao,i,j,'molden')
                        pc.add_tdm(dm1_ao.conj().T,j,i,'molden')
pc.compute_projected_cap()
W=pc.get_projected_cap()

In this example, the files are assumed to located in ./data_files with names cid1trfl.FROMdrt{drtFrom}.state{i}TOdrt{drtTo}.state{i}.iwfmt, which is consistent with them having been generated by the utilities/write_iwfmt.bash script.

Alternatively, one can absolute paths:

dm1_ao = parser.sdm_ao(1,filename='data_files/cid1trfl.FROMdrt1.state1TOdrt1.state1.iwfmt')
pc.add_tdm(dm1_ao,0,0,'molden')

Step 4: Generate and analyze eigenvalue trajectories

H0 and W, or the CAP object can be used to construct a CAPHamiltonian object.

from pyopencap.analysis import CAPHamiltonian
CAPH = CAPHamiltonian(H0=H0,W=W_mat)
# equivalently
CAPH = CAPHamiltonian(pc=pc)

See the analysis section for more details.

Officially supported methods

MR-CISD has been officially tested, though the interface should work with other methods. Please contact us if you find success or have issues using any other methods so we can add official support!