from os import environ
import numpy as np
environ["NUMEXPR_NUM_THREADS"] = "1"
environ["OMP_NUM_THREADS"] = "1"
environ["MKL_NUM_THREADS"] = "1"
[docs]def wf_type(aom):
"""
Checks the topology of the AOMs to obtain the type of wavefunction.
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:returns: A string with the type of wave function ('rest', 'unrest' or 'no').
:rtype: str
"""
if isinstance(aom[0][0][0], float):
return "rest"
elif isinstance(aom[1][0][0], float):
return "no"
elif isinstance(aom[1][0][0][0], float):
return "unrest"
else:
raise NameError("Could not find the type of wave function from the AOMs")
[docs]def find_multiplicity(aom):
"""
Checks the topology of the AOMs to get the difference between alpha and beta electrons.
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:returns: A string with the multiplicity of the calculations.
:rtype: str
"""
if len(aom[0][0]) == len(aom[1][0]):
return "singlet"
elif len(aom[0][0]) == (len(aom[1][0]) + 1):
return "doublet"
elif len(aom[0][0]) == (len(aom[1][0]) + 2):
return "triplet"
elif len(aom[0][0]) == (len(aom[1][0]) + 3):
return "quadruplet"
elif len(aom[0][0]) == (len(aom[1][0]) + 4):
return "quintuplet"
else:
return None
[docs]def load_file(file):
"""
Loads a variable from a file.
:param file: Contains the name or the path containing the variable.
:type file: str
:returns: The file loaded.
:rtype: object
"""
from pickle import load
with open(file, "rb") as f:
file = load(f)
return file
[docs]def save_file(file, save):
"""
Saves a variable to a file.
:param file: The variable to be stored.
:type file: object
:param save: The name or the path where the variable will be saved.
:type save: str
:returns: None
"""
from pickle import dump
with open(save, "wb") as f:
dump(file, f)
[docs]def find_distances(arr, geom):
"""
Collects the distance between the atoms in ring connectivity.
:param arr: Indices of the atoms in ring connectivity.
:type arr: list of int
:param geom: Geometry of the molecule as in mol.atom_coords().
:type geom: numpy.ndarray
:returns: List containing the distances of the members of the ring in Bohrs.
:rtype: list of float
"""
distances = []
for i in range(len(arr)):
coord1 = geom[arr[i] - 1]
coord2 = geom[arr[(i + 1) % len(arr)] - 1]
distances.append(np.linalg.norm(coord1 - coord2) * 0.529177249)
return distances
[docs]def find_dis(arr, aom):
"""
Collects the DIs between the atoms in ring connectivity.
:param arr: Indices of the atoms in ring connectivity.
:type arr: list of int
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:returns: List containing the DIs of the members of the ring.
:rtype: list of float
"""
return [4 * np.trace(np.dot(aom[arr[i] - 1], aom[arr[(i + 1) % len(arr)] - 1])) for i in range(len(arr))]
[docs]def find_di(aom, i, j):
"""
Collects the DI between two atoms.
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:param i: Index of the first atom.
:type i: int
:param j: Index of the second atom.
:type j: int
:returns: DI between the atoms i and j.
:rtype: float
"""
return 2 * np.trace(np.dot(aom[i - 1], aom[j - 1]))
[docs]def find_di_no(aom, i, j):
"""
Collects the DI between two atoms for correlated wavefunctions.
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:param i: Index of the first atom.
:type i: int
:param j: Index of the second atom.
:type j: int
:returns: DI between the atoms i and j.
:rtype: float
"""
return np.trace(np.linalg.multi_dot((aom[1], aom[0][i - 1], aom[1], aom[0][j - 1])))
[docs]def find_lis(arr, aom):
"""
Collects the LIs between the atoms in ring connectivity.
:param arr: Indices of the atoms in ring connectivity.
:type arr: list of int
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:returns: List containing the DIs of the members of the ring.
:rtype: list of float
"""
return [2 * np.trace(np.dot(aom[arr[i] - 1], aom[arr[i] - 1])) for i in range(len(arr))]
[docs]def find_ns(arr, aom):
"""
Collects the atomic populations of all the atoms in the ring.
:param arr: Indices of the atoms in ring connectivity.
:type arr: list of int
:param aom: The Atomic Overlap Matrices (AOMs) in the MO basis.
:type aom: list of matrices
:returns: List containing the atomic populations of the members of the ring.
:rtype: list of float
"""
return [2 * np.trace(aom[arr[i] - 1]) for i in range(len(arr))]
[docs]def av1245_pairs(arr):
"""
Collects the series of atoms that fulfill the 1-2-4-5 relationship for the AV1245 calculation.
:param arr: Indices of the atoms in ring connectivity.
:type arr: list of int
:returns: List containing the series of atoms that fulfill the 1-2-4-5 relationship.
:rtype: list of tuple
"""
return [(arr[i % len(arr)], arr[(i + 1) % len(arr)], arr[(i + 3) % len(arr)], arr[(i + 4) % len(arr)])
for i in range(len(arr))]
[docs]def mol_info(mol=None, mf=None, save=None, partition=None):
"""
Obtains information from the molecule and the calculation to complement the main code function without requiring the 'mol' and 'mf' objects.
:param mol: PySCF Mole object.
:type mol: pyscf.gto.Mole
:param mf: PySCF SCF object.
:type mf: pyscf.scf.hf.SCF
:param save: String with the name of the file to save the information.
:type save: str
:param partition: String with the name of the partition.
:type partition: str
:returns: Dictionary with the information of the molecule and the calculation.
:rtype: dict
"""
info = {}
info.update({"partition": partition})
if mol:
info.update({
"symbols": [mol.atom_symbol(i) for i in range(mol.natm)],
"atom_numbers": [i + 1 for i in range(mol.natm)],
"basisset": mol.basis,
"geom": mol.atom_coords()
})
if mf:
info.update({
"calctype": mf.__class__.__name__,
"energy": mf.e_tot,
"method": mf.__module__
})
if "dft" in mf.__module__ and mf.xc is not None:
info.update({"xc": mf.xc})
else:
info.update({"xc": "None"})
if save:
save_file(info, save)
return info
[docs]def mapping(arr, perm):
"""
Maps the elements of a list according to a permutation.
:param arr: List of elements.
:type arr: list
:param perm: Permutation of the elements.
:type perm: list of int
:returns: List of elements corresponding to a given permutation.
:rtype: list
"""
return [arr[i] for i in range(len(perm))]
[docs]def get_natorbs(mf, S):
"""
Obtains the natural orbitals from the SCF calculation.
:param mf: PySCF SCF object.
:type mf: pyscf.scf.hf.SCF
:param S: Overlap matrix in an AO basis.
:type S: numpy.ndarray
:returns: List containing the occupancies and the Natural Orbitals.
:rtype: tuple of (numpy.ndarray, numpy.ndarray)
"""
from scipy.linalg import eigh
import numpy as np
rdm1 = mf.make_rdm1(ao_repr=True) # In AO basis
# Only Restricted orbitals are supported with current version
if np.ndim(rdm1) == 3:
D = rdm1[0] + rdm1[1]
elif np.ndim(rdm1) == 2:
D = rdm1
else:
raise ValueError(" | Could not find dimensions for the 1-RDM")
occ, coeff = eigh(np.linalg.multi_dot((S, D, S)), b=S)
coeff = coeff[:, ::-1] # Order coefficients
occ = occ[::-1] # Order occupancies
occ[occ < 10 ** -12] = 0.0 # Set small occupancies to 0
occ = np.diag(occ)
return occ, coeff
[docs]def build_eta(mol):
"""
Builds the eta matrices for the partitioning. They consist of a block-truncated matrix with all elements being
zero except for the diagonal elements corresponding to the basis functions of a given atom, which are set to one.
:param mol: PySCF Mole object.
:type mol: pyscf.gto.Mole
:returns: List containing the eta matrices for each atom.
:rtype: list of numpy.ndarray
"""
import numpy as np
eta = [np.zeros((mol.nao, mol.nao)) for _ in range(mol.natm)]
for i in range(mol.natm):
start = mol.aoslice_by_atom()[i, -2]
end = mol.aoslice_by_atom()[i, -1]
eta[i][start:end, start:end] = np.eye(end - start)
return eta
