The program has been written by using parts of the program APOST by I. Mayer and A. Hamza, Budapest, 2000-2003
The numerical integration utilizes the subroutines for Lebedev quadrature downloaded from here. The appropriate reference is: V.I. Lebedev, and D.N. Laikov "A quadrature formula for the sphere of the 131st algebraic order of accuracy" Doklady Mathematics, 59 477-481 1999
The program makes use of Libxc library when necessary, using the F90 interfaces provided by the authors (see http://www.tddft.org/programs/libxc)
We are extremely grateful for the possibility of using these routines!
We also acknowledge R. Oswald for the technical support in code parallelization and in compilation setup preparation
- General structure
- Specific keywords and description
- Input examples
- Extracting .fchk, .dm1 and .dm2 files from pySCF
The code requires the wavefunction information in a Gaussian-type formatted checkpoint file name-input.fchk
and an additional input file name-input.inp requesting the type of calculation.
The formatted checkpoint file can be directly obtained from a Gaussian or Q-Chem run, or can also be generated by I/O converters such as [IOdata] (https://github.com/theochem/iodata).
In the utils folder, we also provide a Python file (apost3d.py) that gathers some interface functions for PySCF, including pyscf2fchk to generate an fchk file from a PySCF run.
The input file is a free format ASCII file that contains one or more input blocks enclosed between a header (# HEADER) and a "#" sign
that gather the desired options for the calculation. The main and mandatory input block is # METHOD, in which both the Atom in Molecule (AIM) and
the type of analysis (one or more) to be performed must be indicated. Additional options like creating cube files or defining fragments are also specified in this block section.
Some WF tools or options require additional block sections for tool-specific options. The order of keywords within the block section is indifferent, but notice they are case-sensitive.
| Keyword | Description |
|---|---|
| MULLIKEN | Hilbert-space Mulliken |
| LOWDIN | Hilbert-space Lowdin |
| LOWDIN-DAVIDSON | Hilbert-space Lowdin-Davidson |
| NAO-BASIS | Hilbert-space based on Natural Atomic Orbitals |
| HIRSH | Real-space Hirshfeld |
| HIRSH-IT | Real-space Hirshfeld iterative |
| BECKE-RHO | Real-space Becke-rho (deprecated) |
| TFVC | Real-space Topological Fuzzy Voronoi Cells |
| Keyword | Description |
|---|---|
| EFFAO | Effective Atomic/Fragment Orbitals from the electron density |
| UEFFAO | Spin-resolved effective Atomic/Fragment Orbitals from the alpha and beta densities |
| EFFAO-U | Effective Atomic/Fragment Orbitals from the paired and unpaired densities |
| EOS | Effective Oxidation States (EOS) analysis |
| EOS-U | Effective Oxidation States analysis from the paired and unpaired density functions (EOS-U) |
| OS-CENTROID | Oxidation states from centroids of localized orbitals |
| OSLO | Oxidation States Localized Orbitals (OSLO). Requires additional # OSLO block section |
| SPIN | Local Spin Analysis (LSA) |
| ENPART | Real-space-only molecular energy decomposition (IQA). Requires additional # ENPART block section |
| EDAIQA | Real-space-only molecular energy decomposition of Energy Decomposition Analysis (EDA) terms. Requires additional # EDAIQA block section |
| POLAR | Bader-Keith decomposition of the dipole moment |
| Keyword | Description |
|---|---|
| DOFRAGS | Definition of molecular fragments for the calculations. Requires additional # FRAGMENTS block section |
| DOINT | Generate *.int files for each atom with the Atomic Overlap Matrices in MO basis for the given AIM. These can be read with ESI program |
| CUBE | Plots cube-type files of the Effective Atomic/Fragment Orbitals. Requires additional # CUBE block section |
| DENS=val | Integer val controls which of the P-matrices present in the fchk file is to be used. Default val=1 |
| QCHEM | Required if the .fchk file originates from a Q-Chem calculation (different ordering of sections within) |
| DM=val | Integer val indicates files with the matrix representation of the RDM1 and RDM2 in MO basis will be provided (only for correlated WF methods). Requires additional # DM block section. If val=1 the RDM1 file will be provided. val=2 indicates both RDM1 and RMD2 files will be provided. These files can be generated using an auxiliary function provided in apost3d.py |
List of supported AIM definitions:
| Keyword | Description |
|---|---|
| MULLIKEN | Hilbert-space Mulliken |
| LOWDIN | Hilbert-space Lowdin |
| LOWDIN-DAVIDSON | Hilbert-space Lowdin-Davidson |
| NAO-BASIS | Hilbert-space based on Natural Atomic Orbitals |
| TFVC | Default AIM |
Extra options:
| Keyword | Description |
|---|---|
| FOLI TOLERANCE=val | Tolerance value for the OSLO iterative procedure. Default val=3 (threshold used = 10^(-val)) |
| BRANCH ITERATION=val | Iteration in which the user invokes a branching. Default val=0 (no branching). This part of the code is currently unavailable |
| PRINT NON-ORTHO | Collect the resulting OSLOs (non-orthogonal) in a new .fchk file. By default, only the final set of orthogonalized OSLOs is provided in a new .fchk file |
Important: TFVC keyword must be included in # METHOD block section, independently of the AIM requested for the OSLO calculation.
List of supported specific electronic structure methods :
| Keyword | Description |
|---|---|
| HF | Decomposition of the Hartree-Fock energy |
| LDA | Decomposition of the LDA energy |
| BP86 | Decomposition of the BP86 energy |
| B3LYP | Decomposition of the B3LYP energy |
| LIBRARY | Decomposition of the molecular energy of a general DFT functional (if supported by Libxc) |
| CASSCF | Decomposition of the CASSCF energy. Requires of the .dm1 and .dm2 files introduced in an additional # DM block section |
| CISD | Decomposition of the CISD energy. Requires of the .dm1 and .dm2 files introduced in an additional # DM block section |
When selecting LIBRARY, the exchange and correlation (or exchange-correlation) functional's ID given by the Libxc library must be provided. To do so, the following additional keywords are required:
| Keyword | Description |
|---|---|
| EXC_FUNCTIONAL=val | Exchange-correlation functional ID |
| EX_FUNCTIONAL=val | Exchange functional ID |
| EC_FUNCTIONAL=val | Correlation functional ID |
Extra options:
| Keyword | Description |
|---|---|
| CORRELATION | For correlated wavefunctions, decomposition of the exchange and correlation energies separately. By default, the code decomposes exchange-correlation altogether |
| THREBOD=val | Integer val sets a threshold for computing exactly the diatomic exchange-correlation terms. If the bond order between a pair of atoms is smaller than the threshold, an approximate multipolar approach is used instead. Default val=100 (actual threshold used val/10000.0d0) |
| ANALYTIC | Perform semianalytical integration of the two-electron energy (decomposition into one-center terms only) (under development) |
| TWOELTOLER=val | Real val sets the threshold (in kcal/mol) for activating the zero-error scheme on the two-electron energy. Default val=0.25d0 |
| MOD-GRIDTWOEL | User-defined integration setup for the two-electron energy. Requires additional # GRID block section |
Important additional information: The zero-error scheme requires separate kinetic, electron-nuclear, and two-electron energies. They can be read from the Gaussian log file (if keyword #p was used), and incorporated in the fchk file using the following utilities before the APOST-3D run as
- For G09 -> $APOST3D_PATH/utils/get_energy mol.log >> mol.fchk
- For G16 -> $APOST3D_PATH/utils/get_energy_g16 mol.log >> mol.fchk
For .fchk files obtained from pySCF no action is required.
| Keyword | Description |
|---|---|
| MAX_OCC=val | Maximal EFO occupation for generating its cube files. Integer val is in the [0, 1000] range and the occupation is defined as val/1000.0d0) |
| MIN_OCC=val | Minimal EFO occupation for generating the cube file. val is an integer in the [0, 1000] range and the occupation is defined as val/1000.0d0) |
| Keyword | Description |
|---|---|
| RADIAL=val | Integer number of radial grid points. Default val=150 |
| ANGULAR=val | Integer number of angular grid points according to the Lebedev-Laikov spherical grids. Default val=590. Other recommended values: 146, 434, or 974 |
| rr00=val | Double-precision variable that sets the distance that contains half of the radial integration points. It corresponds to the value of rm in Eq. 25 of [Becke's multicenter integration scheme] (https://doi.org/10.1063/1.454033). Default val=0.5d0 |
| phb1=val | Double-precision value for the first rotation angle (in radians) of the grid for the second electron coordinates. Default val=0.169 |
| phb2=val | Double-precision value for the second rotation angle (in radians) for the zero-error strategy). Default val=0.170 |
Important additional information: The default values of phb1 and phb2 have been optimized for a (150,590) grid. Changing the radial or angular points would require recalibration of the phb1 and phb2 values.
Include in the first line the filename of the RDM1 (.dm1). If DM=2 was set in METHOD block section, include in the second line the filename of the RDM2 (.dm2). If the dm1 and dm2 files were generated by apost3d.py utility, include in a new line the keyword pySCF
Definition of fragments using their ordering in the fchk file. For instance
# FRAGMENTS
3 ## Number of fragments ##
1 ## Number of atoms for fragment 1 ##
1 ## Atom list for fragment 1 ##
2 ## Number of atoms for fragment 2 ##
3 5 ## Atom list for fragment 2 ##
3 ## Number of atoms for fragment 3 ##
2 4 6 ## Atom list for fragment 3 ##
#
In this case, 3 fragments have been defined. the first consists solely of one center (atom 1), the second is formed by two atoms (3 and 5) and the third is formed by three centers (2,4 and 6). Notice that all atoms in the fchk file must be assigned to a fragment. All remaining atoms can be assigned to the last fragment by writing -1 instead of its number of atoms.
Sample input files are provided with the most common types of analyses currently implemented in APOST-3D. Notice that more than one analysis tool using the AIM can be requested with the same input file.
# METHOD
NAO-BASIS
EOS
DOFRAGS
CUBE
#
# FRAGMENTS
2
1
1
-1
#
# CUBE
MAX_OCC=700
MIN_OCC=300
#
In this example, the evaluation of the effective oxidation states is requested, together with the cube generation of the effective fragment orbitals with net occupancies within the [0.300,0.700] range.
The AIM requested is NAO-BASIS. which requires of an additional inputname.nao file with the transformation matrix from the original AO to NAO basis. This file corresponds to the FILE.33 file generated, for instance, with Gaussian coupled to the NBO software. For this aim, in a Gaussian input file the user must add the pop=(full,nboread) keyword, together with the $NBO AONAO=W $END additional line at the end of the file. This generates the FILE.33 file to be renamed as inputname.nao.
Fragments have been defined. In particular, the system is split into two fragments: fragment 1 consists of atom 1, and fragment 2 gathers the rest of atoms of the molecule.
# METHOD
TFVC
ENPART
DOFRAGS
#
# ENPART
LIBRARY
EX_FUNCTIONAL=106
EC_FUNCTIONAL=132
THREBOD=-1
MOD-GRIDTWOEL
#
# FRAGMENTS
2
1
1
-1
#
# GRID
RADIAL 150
ANGULAR 590
rr00 0.5
phb1 0.169
phb2 0.170
#
In this example, the IQA decomposition of the KS-DFT energy is requested. Functional ID=106 (B88) and ID=132 (P86) correspond to the well-known BP86 functional.
The AIM requested is TFCV and fragments have been defined. In particular, the system is split into two fragments: fragment 1 consists of atom 1, and fragment 2 gathers the rest of atoms of the molecule. The diatomic xc energy components are computed exactly for all atomic pairs. User-specific options for the decomposition of the two-electron energy are given in GRID block section (in this case they coincide with the default strongly recommended values).
# METHOD
TFVC
SPIN
DM 2
#
# DM
filename.dm1
filename.dm2
pySCF
#
In this example, the local spin analysis of a correlated WF (e.g. CASSCF) from a pySCF run is requested.
The files containing the information of the RDM1 and RDM2 (.dm1 and .dm2) are mandatory for this type of analysis. They can be obtained from a pySCF run using the apost3d.py utility, as described here.
The DM section contains the name of the two files (with extension), together with the keyword of the code they have been created with.
The AIM requested is TFVC and fragments have not been defined.
# METODE
TFVC
OSLO
DOFRAGS
#
# OSLO
LOWDIN
FOLI TOLERANCE 3
BRANCH ITERATION 0
PRINT NON-ORTHO
#
# FRAGMENTS
3
1
1
2
2 3
-1
#
In this example, the evaluation of the oxidation states using the OSLO method is requested.
The AIM requested for the numerical integration is TFVC (mandatory for OSLO calculation) and fragments have been defined. In particular, the system is split into three fragments: fragment 1 consists of atom 1, fragment 2 consists of two atoms (2 and 3), and fragment 3 gathers the rest of atoms of the molecule.
The AIM requested for FOLI values evaluation is LOWDIN, and the tolerance value for its selection is 3 (real(10^-3), default value). Printing of the selected OSLOs pre-orthogonalization (non-orthogonal if two or more OSLOs are selected in the same iteration AND belong to different fragments) is invoked. As a result, a new .fchk file is created containing these OSLOs.
# METODE
TFVC
EOS
SPIN
DOFRAGS
#
# FRAGMENTS
2
1
1
4
3 4 5 6
#
In this example, the evaluation of the effective oxidation states and the local spin analysis. For correlated (multireference) WFs, the local spin analysis requires the information of the RDM1 and RDM2 (see Example 3).
The AIM requested is TFVC and fragments have been defined. In particular, the system is split into two fragments: fragment 1 consists of atom 1, and fragment 2 of atoms 2-6 (the system has 6 atoms).
A battery of real-case examples of calculations are provided in the $APOST3D_PATH/input_examples folder, together with the output files generated by execution of the code.
The apost3d.py file provided in $APOST3D_PATH/utils folder collects several functions for creating .fchk, .dm1 and .dm2 files from a pySCF run. The location of the file must be added to the PYTHONPATH variable, for instance by doing
export PYTHONPATH=$PYTHONPATH:$APOST3D_PATH/utils
As example, we provide an input of a CASSCF(2,2) calculation for the HF molecular system in the singlet spin state.
from pyscf import gto, scf, mcscf
import apost3d as apost
molname = 'HF-CASSCF'
mol=gto.M()
mol.atom='''
H 0.0000000 0.0000000 0.0000000
F 0.0000000 0.0000000 0.9500000
'''
mol.basis='def2tzvp'
mol.spin = 0
mol.charge = 0
mol.cart= False
mol.symmetry = False
mol.verbose=4
mol.max_memory=40000
mol.build()
# Calculate RHF reference (guess) #
mf = scf.RHF(mol)
mf.kernel()
# Once our guess (converged HF) is created, time to create and converge the CASSCF #
mycas = mcscf.CASSCF(mf,2,2) # Active space size: Orbitals, Electrons #
cas_list = [5,6] # Pick orbitals for CAS space, 1-based indices #
mo = mcscf.sort_mo(mycas, mf.mo_coeff, cas_list) # Just orders the MOs #
mycas.natorb = True
# Run CAS #
mycas.kernel(mo)[0]
# Resulting orbitals printing in .fchk file #
s = mf.get_ovlp()
apost.write_fchk(mol, mycas, molname, s)
# .dm1 and .dm2 files creation ONLY for multiconfig WFs #
apost.write_dm12(mol, mycas, molname)
For single-determinant calculations (Hartree-Fock or KS-DFT), the user only requires the .fchk file for running APOST-3D. The write_fchk function is independent of the WF-type, being thus called as in the example input provided for CASSCF.
For illustrative purposes, one can create the .fchk file from a Hartree-Fock calculation as
from pyscf import gto, scf
import apost3d as apost
molname = 'HF-Hartree-Fock'
mol=gto.M()
mol.atom='''
H 0.0000000 0.0000000 0.0000000
F 0.0000000 0.0000000 0.9500000
'''
mol.basis='def2tzvp'
mol.spin = 0
mol.charge = 0
mol.cart= False
mol.symmetry = False
mol.verbose=4
mol.max_memory=40000
mol.build()
# Calculate RHF #
mf = scf.RHF(mol)
mf.kernel()
# Resulting orbitals printing in .fchk file #
s = mf.get_ovlp()
apost.write_fchk(mol, mycas, molname, s)