chemtools.basisset module

class BasisSet(name, element, family=None, kind=None, functions=None, info=None)[source]

Bases: object

Basis set module supporting basic operation on basis sets.

Args:
name : str
Name of the basis set
element : str
Symbol of the element
Kwargs:
kind : str
Classification of the basis set functions, diffuse, tight
family : str
basis set family
functions : dict
Dict of functions with s, p, d, f, ... as keys
params : list of dicts
Parameters for generating the functions according to the model
append(other)[source]

Append functions from another BasisSet object

Args:
other : BasisSet
BasisSet object whose functions will be added to the existing ones
completeness_profile(zetas)[source]

Calculate the completeness profile of each shell of the basis set

Args:
zetas : numpy.array
Scaning exponents
Returns:
out : numpy.array
numpy array with values of the profile (shells, zetas)
contraction_matrix(shell)[source]

Return the contraction coefficients for a given shell in a matrix form with size ne * nc, where ne is the number of exponents and nc is the number of contracted functions

Args:
shell : str
shell label, s, p, d, ...
Returns:
out : 2D numpy.array
2D array with contraction coefficients
contraction_scheme()[source]

Return a string describing the contraction scheme.

contraction_type()[source]

Try to determine the contraction type: segmented, general, uncontracted, unknown.

contractions_per_shell()[source]

Calculate how many contracted functions are in each shell.

Returns:
out : list of ints
classmethod from_file(fname=None, fmt=None, name=None)[source]

Read and parse a basis set from file and return a BasisSet object

Args:
fname : str
File name
fmt : str
Format of the basis set in the file (molpro, gamessus)
name : str
Name of the basis set
Returns:
out : BasisSet or dict
Basisset object parsed from file or dictionary of BasisSet objects
classmethod from_json(jsonstring, **kwargs)[source]

Instantiate the BasisSet object from a JSON string

Args:
jsonstring: str
A JSON serialized string with the basis set
classmethod from_optpars(x0, funs=None, name=None, element=None, explogs=False)[source]

Return a basis set object generated from a sequence based on the specified arguments.

Args:
x0 : list or numpy.array
Parameters to generate the basis set as a continuous list or array
funs : list of tuples
A list of tuple specifying the shell type, number of functions and parameters, e.g. [(‘s’, ‘et’, 4, (0.5, 2.0)), (‘p’, ‘et’, 3, (1.0, 3.0))]
name : str
Name of the basis set
element : str
Chemical symbol of the element
Returns:
out : BasisSet
static from_pickle(fname, **kwargs)[source]

Read a pickled BasisSet object from a pickle file

Args:
fname : str
File name containing the BasisSet
kwargs : dict
Extra arguments for the pickle.load method
Raises:
UnicodeDecodeError
When you try to read a python2 pickle using python3, to fix that use encoding=’latin1’ option
classmethod from_sequence(funs=None, name=None, element=None)[source]

Return a basis set object generated from a sequence based on the specified arguments.

Args:
funs : list of tuples
A list of tuple specifying the shell type, number of functions and parameters, e.g. [(‘s’, ‘et’, 4, (0.5, 2.0)), (‘p’, ‘et’, 3, (1.0, 3.0))]
name : str
Name of the basis set
element : str
Chemical symbol of the element
Returns:
out : BasisSet
classmethod from_str(string, fmt=None, name=None)[source]

Parse a basis set from string

Args:
string : str
A string with the basis set
fmt : str
Format of the basis set in the file: molpro, gamessus
name : str
Name of the basis set
Returns:
out : BasisSet or dict
Basisset object parsed from string or dictionary of BasisSet objects
get_exponents(asdict=True)[source]

Return the exponents of a given shell or if the shell isn’t specified return all of the available exponents

Args:
asdict (bool): if True a dict with exps per shell is
returned
get_filename(ext=None)[source]
nf(spherical=True)[source]

Calculate the number of basis functions

Args:
spherical : bool
flag indicating if spherical or cartesian functions should be used, default: True
Returns:
res : int
number of basis functions
normalization()[source]

For each function (contracted) calculate the norm and return a list of tuples containing the shell, function index and the norm respectively.

normalize()[source]

Normalize contraction coefficients for each contracted functions based on the primitive overlaps so that the norm is equal to 1.

nprimitive(spherical=True)[source]

Return the number of primitive functions assuming sphrical or cartesian gaussians.

Args:
spherical : bool
A flag to select either spherical or cartesian gaussians
Returns:
out : int
Number of primitive function in the basis set
partial_wave_expand()[source]

From a given basis set with shells spdf... return a list of basis sets that are subsets of the entered basis set with increasing angular momentum functions included [s, sp, spd, spdf, ...]

primitives_per_contraction()[source]

Calculate how many primities are used in each contracted function.

Returns:
out : list of ints
primitives_per_shell()[source]

Calculate how many primitive functions are in each shell.

Returns:
out : list of ints
print_functions(efmt='20.10f', cfmt='15.8f')[source]

Return a string with the basis set.

Args:
efmt : str
string describing output format for the exponents, default: “20.10f”
cfmt : str
string describing output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string
shell_overlap(shell)[source]

Calculate the overlap integrals for a given shell

Args:
shell : str
Shell
Returns:
out : numpy.array
Overlap integral matrix
sort(reverse=False)[source]

Sort shells in the order of increasing angular momentum and for each shell sort the exponents.

Args:
reverse : bool
If False sort the exponents in each shell in the descending order (default), else sort exponents in ascending order
to_cfour(comment='', efmt='15.8f', cfmt='15.8f')[source]

Return a string with the basis set in (new) CFOUR format.

Args:
comment : str
comment string
efmt : str
string describing output format for the exponents, default: “20.10f”
cfmt : str
string describing output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string in Cfour format
to_dalton(fmt='prec')[source]

Return a string with the basis set in DALTON format.

Args:
fmt : str
string describing output format for the exponents and coefficents
  • prec “20.10f”
  • default “10.4f”
  • or python format string e.g. “15.8f”
Returns:
res : str
basis set string in Dalton format
to_gamessus(efmt='20.10f', cfmt='15.8f')[source]

Return a string with the basis set in GAMESS(US) format.

Args:
efmt : str
string describing output format for the exponents, default: “20.10f”
cfmt : str
string describing output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string in Gamess(US) format
to_gaussian(efmt='20.10f', cfmt='15.8f')[source]

Return a string with the basis set in Gaussian format.

Args:
efmt : str
string describing output format for the exponents, default: “20.10f”
cfmt : str
string describing output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string in Gaussian format
to_json(fname=None, **kwargs)[source]

Serizalize the basis set object to JSON format

to_latex(efmt='20.10f', cfmt='15.8f')[source]

Return a string with the basis set as LaTeX table/

Args:
efmt : str
Output format for the exponents, default: “20.10f”
cfmt : str
Output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string in LaTeX format
to_molpro(withpars=False, efmt='20.10f', cfmt='15.8f')[source]

Return a string with the basis set in MOLPRO format.

Args:
withpars : bool
A flag to indicate whether to wrap the basis with basis={ } string
efmt : str
Output format for the exponents, default: “20.10f”
cfmt : str
Output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string
to_nwchem(efmt='20.10f', cfmt='15.8f')[source]

Return a string with the basis set in NWChem format.

Args:
efmt : str
string describing output format for the exponents, default: “20.10f”
cfmt : str
string describing output format for the contraction coefficients, default: “15.8f”
Returns:
res : str
basis set string in NwChem format
to_pickle(fname=None)[source]

Save the basis set in pickle format under the filename fname

Args:
fname : str
File name
uncontract(copy=False)[source]

Uncontract the basis set. This replaces the contraction coefficients in the current object.

Args:
copy : bool
If True return an uncontracted copy of the basis set rather than uncontracting in place, default is False.
eventemp(numexp, params)[source]

Generate a sequence of nf even tempered exponents accodring to the even tempered formula

\[\zeta_i = \alpha \cdot \beta^{i-1}\]
Args:
numexp : int
Number fo exponents to generate
params : tuple of floats
Alpha and beta parameters
Returns:
res : numpy array
Array of generated exponents (floats)
generate_exponents(formula, nf, params)[source]

Generate a sequence of exponents from a specified formula

Args:
formula : str
name of the sequence from which the exponents are generated, one of:
  • et, even, eventemp, even tempered
  • wt, well, welltemp, well tempered
  • le, legendre
  • exp, exponents
nf : int
number of exponents
params : tuple
parameters for the formula
Returns:
numpy.array with exponents
get_num_params(ftuple)[source]

Get the number of parameter for a given 4-tuple desciribing functions

>>> get_num_params(('s', 'et', 20, (0.4, 2.0)))
2
has_consecutive_indices(shell)[source]

Check if all the contracted functions have consecutive indices

Args:
shell : dict
Basis functions for a given shell asa dict with structure {'e' : np.array(), 'cf': [np.array(), np.array(), ...]}
legendre(numexp, coeffs)[source]

Generate a sequence of nf exponents from expansion in the orthonormal legendre polynomials as described in: Peterson, G. A. et.al J. Chem. Phys., Vol. 118, No. 3 (2003), eq. (7).

\[\ln \zeta_i = \sum^{k_{\max}}_{k=0} A_k P_k \left(\frac{2j-2}{N-1}-1\right)\]
Args:
numexp : int
Number fo exponents to generate
params : tuple of floats
Polynomial coefficients (expansion parameters)
Returns:
res : numpy array
Array of generated exponents (floats)
merge(first, other)[source]

Merge functions from two BasisSet objects

Args:
first : BasisSet
First BasisSet object to merge
other : BasisSet
Second BasisSet object to merge
Returns:
our : BasisSet
BasisSet instance with functions from first and other merged
ncartesian(l)[source]

Calculate the number of cartesian components of a function with a given angular momentum value l.

nspherical(l)[source]

Calculate the number of spherical components of a function with a given angular momentum value l.

primitive_overlap(l, a, b)[source]

Calculate the overlap integrals for a given shell l and two sets of exponents

\[S(\zeta_i, \zeta_j) = \frac{2^{l + \frac{3}{2}} \left(\zeta_{1} \zeta_{2}\right)^{\frac{l}{2} + \frac{3}{4}}}{\left(\zeta_{1} + \zeta_{2}\right)^{l + \frac{3}{2}}}\]
Args:
l : int
Angular momentum quantum number of the shell
a (M,) : numpy.array
First vector of exponents
b (N,) : numpy.array
Second vector of exponents
Returns:
out (M, N) : numpy.array
Overlap integrals
reorder_shell_to_consecutive(shell)[source]

Reorder the exponents so that the indices of the contracted functions have consecutive inidices.

Args:
shell : dict
Basis functions for a given shell asa dict with structure {'e' : np.array(), 'cf': [np.array(), np.array(), ...]}
Returns:
shell : dict
Same shell as on input but with reordered exponents and relabelled contracted functions
sliceinto(toslice, chunk_sizes)[source]

Slice a list toslice into chunks of sizes defined in chunk_sizes

>>> sliceinto(list(range(10)), (6, 4))
[[0, 1, 2, 3, 4, 5], [6, 7, 8, 9]]
splitlist(l, n)[source]

Split a list into sublists of size n

welltemp(numexp, params)[source]

Generate a sequence of nf well tempered exponents accodring to the well tempered fromula

\[\zeta_i = \alpha \cdot \beta^{i-1} \cdot \left[1 + \gamma \cdot \left(\frac{i}{N}\right)^{\delta}\right]\]
Args:
numexp : int
Number fo exponents to generate
params : tuple of floats
Alpha, beta, gamma and delta parameters
Returns:
res : numpy.array
Array of generated exponents (floats)
xyzlist(l)[source]

Generate an array of \(l_x\), \(l_y\), \(l_z\) components of cartesian gaussian with a given angular momentum value in canonical order.

For exampe:
  • l = 0 generates the result array([[0, 0, 0]])
  • l = 1 generates array([[1, 0, 0], [0, 1, 0], [0, 0, 1])
  • etc.

The functions are coded by triples of powers of x, y, z, namely [1, 2, 3] corresponds to \(xy^{2}z^{3}\).

Args:
l : int
Angular momentum value
Returns:
out ((l+1)*(l+2)/2, 3) : numpy.array
Array of monomial powers, where columns correspond to x, y, and z respectively and rows correspond to functions.
zetas2eventemp(zetas)[source]

From a set of exponents (zetas), using least square fit calculate the approximate $alpha$ and $beta$ parameters from the even tempered expansion.

Args:
zetas : list
list of exponents (floats) to be fitted
Returns:
coeff : numpy.array
numpy array of legendre expansion coeffcients of length kmax
zetas2legendre(zetas, kmax)[source]

From a set of exponents (zetas), using least square fit calculate the expansion coefficients into the legendre polynomials of the order kmax.

Args:
kmax : int
length of the legendre expansion
zetas : list
list of exponents (floats) to be fitted
Returns:
coeff : numpy.array
numpy array of legendre expansion coeffcients of length kmax
zlmtoxyz(l)[source]

Generates the expansion coefficients of the real spherical harmonics in terms of products of cartesian components. Method based on [1]

[1]Schlegel, H. B., & Frisch, M. J. (1995). “Transformation between Cartesian and pure spherical harmonic Gaussians”. International Journal of Quantum Chemistry, 54(2), 83–87. doi:10.1002/qua.560540202
Args:
l : int
Angular momentum value
Returns:
out \(((l+1)*(l+2)/2, 2*l + 1)\) : numpy.array
Expansion coefficients of real spherical harmonics in terms of cartesian gaussians