API Reference¶
Thermochemistry¶
-
class
panther.thermochemistry.
Thermochemistry
(vibenergies, atoms, *args, **kwargs)[source]¶ Calculate thermochemistry in harmonic approximation
Parameters: vibenergies : numpy.array
Vibrational energies in Joules
atoms : ase.Atoms
Atoms obect
phase : str
Phase, should be either gas or solid
pointgroup : str
symmetrynumber : str
If pointgroup is specified symmetrynumber is obsolete, since it will be inferred from the pointgroup
-
get_internal_energy
(T=273.15)[source]¶ Return the internal energy U
Parameters: T : float
Temperature in K
-
get_qvibrational
(T=273.15, uselog=True)[source]¶ Calculate the vibrational partition function at temperature T in kJ/mol
Parameters: T : float
Temperature in K
uselog : bool
When True return the natural logarithm of the partition function
Notes
\[q_{vib}(T) = \prod^{3N-6}_{i=1}\frac{1}{1 - \exp(-h\omega_{i}/k_{B}T)}\]
-
get_vibrational_energy
(T=273.15)[source]¶ Calculate the vibational energy correction at temperature T in kJ/mol
Parameters: T : float
Temperature in K
Notes
\[U_{vib}(T) = \frac{R}{k_{B}}\sum^{3N-6}_{i=1} \frac{h\omega_{i}}{\exp(h\omega_{i}/k_{B}T) - 1}\]
-
get_vibrational_entropy
(T=273.15)[source]¶ Calculate the vibrational entropy at temperature T in kJ/mol
Parameters: T : float
Temperature in K
Notes
\[S_{vib}(T) = R\sum^{3N-6}_{i=1}\left[ \frac{h\omega_{i}}{k_{B}T(\exp(h\omega_{i}/k_{B}T) - 1)} - \ln(1 - \exp(-h\omega_{i}/k_{B}T)) \right]\]
-
get_vibrational_heat_capacity
(T=273.15)[source]¶ Return the heat capacity
Parameters: T : float
Temperature in K
Notes
\[C_{p,vib}(T) = R\sum^{3N-6}_{i=1} \left(\frac{h\omega_{i}}{k_{B}T}\right)^{2} \frac{\exp(-h\omega_{i}/k_{B}T)}{\left[1 - \exp(-h\omega_{i}/k_{B}T)\right]^{2}}\]
-
NormalModeBFGS¶
-
class
panther.nmrelaxation.
NormalModeBFGS
(atoms, phase, hessian=None, hessian_update='BFGS', logfile='-', trajectory=None, restart=None, proj_translations=True, proj_rotations=True, master=None, verbose=False)[source]¶ Normal mode optimizer with approximate hessian update
Parameters: atoms : ase.Atoms
Atoms object with the structure to optimize
phase : str
Phase, should be either gas or solid
hessian : array_like (N, N)
Initial hessian matrix in eV/Angstrom^2
hessian_update : str
Name of the approximate formula to udpate hessian, one of: BFGS, SR1, DFP
proj_translations : bool
If
True
translational degrees of freedom will be projected from the hessianproj_rotations : bool
If
True
rotational degrees of freedom will be projected from the hessianlogfile : str
Name log the log file
trajectory : str
Name of the trajectory file
-
run
(fmax=0.05, steps=100000000)[source]¶ Run structure optimization algorithm.
This method will return when the forces on all individual atoms are less than fmax or when the number of steps exceeds steps.
-
panther.anharmonicity¶
Methods for solving the one dimentional vibrational eigenproblem
-
panther.anharmonicity.
anharmonic_frequencies
(atoms, T, coeffs, modeinfo)[source]¶ Calculate the anharmonic frequencies
Parameters: atoms : ase.Atoms
Atoms object
T : float
Temperature in K
coeffs : pandas.DataFrame
modeinfo : pandas.DataFrame
-
panther.anharmonicity.
factsqrt
(m, n)[source]¶ Return a factorial like constant
Parameters: m : int
Argument of the series
n : int
Length of the series
Notes
\[f(m, n) = \prod^{n - 1}_{i = 0} \sqrt{m - i}\]
-
panther.anharmonicity.
get_anh_state_functions
(eigenvals, T)[source]¶ Calculate the internal energy
U
and entropyS
for an anharmonic vibrational mode with eigenvalueseigvals
at temperatureT
in kJ/mol\[ \begin{align}\begin{aligned}U &= N_{A}\frac{\sum^{n}_{i=1} \epsilon_{i}\exp(\epsilon_{i}/k_{B}T) }{\sum^{n}_{i=1} \exp(\epsilon_{i}/k_{B}T)}\\S &= N_{A}k_{B}\log(\sum^{n}_{i=1} \exp(\epsilon_{i}/k_{B}T)) + \frac{N_{A}}{T}\frac{\sum^{n}_{i=1} \epsilon_{i}\exp(\epsilon_{i}/k_{B}T) }{\sum^{n}_{i=1} \exp(\epsilon_{i}/k_{B}T)}\end{aligned}\end{align} \]Parameters: eigenvals : numpy.array
Eigenvalues of the anharmonic 1D Hamiltonian in Joules
T : float
Temperature in K
Returns: (U, S) : tuple of floats
Tuple with the internal energy and entropy in kJ/mol
-
panther.anharmonicity.
get_hamiltonian
(rank, freq, mass, coeffs)[source]¶ Compose the Hamiltonian matrix for the anharmonic oscillator with the potential described by the sixth order polynomial.
Parameters: rank : int
Rank of the Hamiltonian matrix
freq : float
Fundamental frequency in hartrees
mass : float
Reduced mass of the mode
coeffs : array
A one dimensional array with polynomial coeffients
Notes
\[H_{ij} = \left\langle \Psi_{i} \left| \hat{H} \right| \Psi_{j} \right\rangle\]where
\[\hat{H} = -\frac{\hbar^2}{2}\frac{\partial^2}{\partial \boldsymbol{Q}^2} + \sum_{\mu=0}^{6}c_{\mu}\boldsymbol{Q}^{\mu}\]and \(\Psi_{i}\) are the standard harmonic oscillator functions.
-
panther.anharmonicity.
harmonic_df
(modeinfo, T)[source]¶ Calculate per mode contributions to the thermodynamic functions in the harmonic approximation
Parameters: modeinfo : pandas.DataFrame
T : float
Temperature in K
Returns: df : pandas.DataFrame
-
panther.anharmonicity.
merge_vibs
(anh6, anh4, harmonic, verbose=False)[source]¶ Form a DataFrame with the per mode thermochemical contributions from three separate dataframes with sixth order polynomial fitted potentia, fourth order fitted potential and harmonic frequencies.
Parameters: anh6 : pandas.DataFrame
anh4 : pandas.DataFrame
harmonic : pandas.DataFrame
Returns: df : pandas.DataFrame
panther.displacements module¶
-
panther.displacements.
calculate_displacements
(atoms, hessian, freqs, normal_modes, npoints=4, modes='all')[source]¶ Calculate displacements in internal coordinates
Parameters: atoms : ase.Atoms
Atoms object with the equilibrium structure
hessian : array_like
Hessian matrix in atomic units
freqs : array_like
Frequencies (square roots of the hessian eigenvalues) in atomic units
normal_modes : array_like
Normal modes in atomic units
npoints : int
Number of points to displace structure, the code will calculate
2*npoints
displacements since + and - directions are takenmodes : str or list/tuple of ints, default ‘all’
Range of the modes for which the displacements will be calculated
Returns: images : dict of dicts
A nested (ordred) dictionary with the structures with mode, point as keys, where point is a number from -4, -3, -2, -1, 1, 2, 3, 4
mi : pandas.DataFrame
DataFrame with per mode characteristics, displacements, masses and vibrational population analysis
-
panther.displacements.
get_internals_and_bmatrix
(atoms)[source]¶ internals is a numpy record array with ‘type’ and ‘value’ records bmatrix is a numpy array n_int x n_cart
Parameters: atoms : ase.Atoms
Atoms object
-
panther.displacements.
get_modeinfo
(hessian, freqs, ndof, Bmatrix_inv, Dmatrix, mwevecs, npoints, internals)[source]¶ Compose a DataFrame with information about the vibrations, each mode corresponds to a separate row
-
panther.displacements.
get_nvibdof
(atoms, proj_rotations, proj_translations, phase, include_constr=False)[source]¶ Calculate the number of vibrational degrees of freedom
Parameters: atoms : ase.Atoms
proj_translations : bool
If
True
translational degrees of freedom will be projected from the hessianproj_rotations : bool
If
True
rotational degrees of freedom will be projected from the hessianinclude_constr : bool
If
True
the constraints will be includedReturns: nvibdof : float
Number of vibrational degrees of freedom
-
panther.displacements.
vib_population
(hessian, freqs, Bmatrix_inv, Dmatrix, internals, mi)[source]¶ Calculate the vibrational population analysis
Parameters: hessian : array_like
Hessian matrix
freqs : array_like
A vector of frequencies (square roots fof hessian eigenvalues)
Bmatrix_inv : array_like
Inverse of the B matrix
Dmatrix : array_like
D matrix
internals : array_like
Structured array with internal coordinates
mi : pandas.DataFrame
Modeinfo
Returns: mi : pandas.DataFrame
Modeinfo DataFrame updated with columns with vibrational population analysis results
panther.io module¶
Module providing functions for reading the input and other related files
-
panther.io.
get_symmetry_number
(pointgroup)[source]¶ Return the symmetry number for a given point group
See also
C. J. Cramer, Essentials of Computational Chemistry, Theories and Models, 2nd Edition, p. 363
Parameters: pointgroup : str
Symbol of the point group
-
panther.io.
parse_arguments
()[source]¶ Parse the input/config file name from the command line, parse the config and return the parameters.
-
panther.io.
print_mode_thermo
(df, info=False)[source]¶ After calculating all the anharmonic modes print the per mode themochemical functions
-
panther.io.
print_modeinfo
(mi, output=None)[source]¶ Print the vibrational population data
Parameters: mi : pandas.DataFrame
output : str
Name of the file to store the printout, if
None
stdout will be used
-
panther.io.
read_bmatdat
()[source]¶ Read the
bmat.dat
file with internal coordiantes and the B matrix produced by the originalwriteBmat
codeReturns: internals, Bmatrix : tuple
Internal coordiantes and B matrix
-
panther.io.
read_em_freq
(fname)[source]¶ Read the file
fname
with the frequencies, reduced masses and fitted fitted coefficients for the potential into a pandas DataFrame.Parameters: fname : str
Name of the file with PES
-
panther.io.
read_pes
(fname)[source]¶ Parse the file with the potential energy surface (PES) into a dict of numpy arrays with mode numbers as keys
Parameters: fname : str
Name of the file with PES
-
panther.io.
read_poscars
(filename)[source]¶ Read POSCARs file with the displaced structures and return an OrderedDict with the Atoms objects
-
panther.io.
read_vasp_hessian
(outcar='OUTCAR', symmetrize=True, convert_to_au=True, dof_labels=False)[source]¶ Parse the hessian from the VASP
OUTCAR
file into a numpy arrayParameters: outcar : str
Name of the VASP output, default is
OUTCAR
symmetrize : bool
If
True
the hessian will be symmetrizedconvert_to_au : bool
If
True
convert the hessian to atomic units, in the other case hessian is returned in [eV/Angstrom**2]dof_labels : bool, default is False
If
True
a list of labels corresponding to the degrees of freedom will also be returnedReturns: hessian : numpy.array
Hessian matrix
Notes
Note
By default VASP prints negative hessian so the elements should be multiplied by -1 to restore the original hessian, this is done by default, hessian in the XML file is NOT symmetrized by default
-
panther.io.
read_vasp_hessian_xml
(xml='vasprun.xml', convert_to_au=True, stripmass=True)[source]¶ Parse the hessian from the VASP
vasprun.xml
file into a numpy arrayParameters: xml : str
Name of the VASP output, default is
vasprun.xml
convert_to_au : bool
If
True
convert the hessian to atomic units, in the other case hessian is returned in [eV/Angstrom**2]dof_labels : bool, default is False
If
True
a list of labels corresponding to the degrees of freedom will also be returnedstripmass : bool
If
True
use VASP default masses to transform hessian to non-mass-weighted formReturns: hessian : numpy.array
Hessian matrix
Notes
Note
By default VASP prints negative hessian so the elements should be multiplied by -1 to restore the original hessian, this is done by default, hessian in the XML file is symmetrized by default
panther.nmrelaxation module¶
-
class
panther.nmrelaxation.
NormalModeBFGS
(atoms, phase, hessian=None, hessian_update='BFGS', logfile='-', trajectory=None, restart=None, proj_translations=True, proj_rotations=True, master=None, verbose=False)[source] Bases:
ase.optimize.optimize.Optimizer
,object
Normal mode optimizer with approximate hessian update
Parameters: atoms : ase.Atoms
Atoms object with the structure to optimize
phase : str
Phase, should be either gas or solid
hessian : array_like (N, N)
Initial hessian matrix in eV/Angstrom^2
hessian_update : str
Name of the approximate formula to udpate hessian, one of: BFGS, SR1, DFP
proj_translations : bool
If
True
translational degrees of freedom will be projected from the hessianproj_rotations : bool
If
True
rotational degrees of freedom will be projected from the hessianlogfile : str
Name log the log file
trajectory : str
Name of the trajectory file
-
log
(grad, grad_nm, step_nm)[source] Print a line with convergence information
-
log_header
()[source] Header for the log with convergence information
-
run
(fmax=0.05, steps=100000000)[source] Run structure optimization algorithm.
This method will return when the forces on all individual atoms are less than fmax or when the number of steps exceeds steps.
-
step
(grad)[source] Calculate the step in cartesian coordinates based on the step in normal modes in the rational function approximation (RFO)
- Args:
- grad : array_like (N,)
- Current gradient
-
update_hessian
(coords, grad)[source] Perform hessian update
Parameters: coords : array_like (N,)
Current coordinates as vector
grad : array_like (N,)
Current gradient
-
-
panther.nmrelaxation.
nmoptimize
(atoms, hessian, calc, phase, proj_translations=True, proj_rotations=True, gtol=1e-05, verbose=False, hessian_update='BFGS', steps=100000)[source]¶ Relax the strcture using normal mode displacements
Parameters: atoms : ase.Atoms
Atoms object with the structure to optimize
hessian : array_like
Hessian matrix in eV/Angstrom^2
calc : ase.Calculator
ASE Calcualtor instance to be used to calculate forces
phase : str
Phase, ‘solid’ or ‘gas’
gtol : float, default=1.0e-5
Energy gradient threshold
hessian_update : str
Approximate formula to update hessian, possible values are ‘BFGS’, ‘SR1’ and ‘DFP’
steps : int
Maximal number of iteration to be performed
verbose : bool
If
True
additional debug information will be printedNotes
Internally eV and Angstroms are used.
See also
Bour, P., & Keiderling, T. A. (2002). Partial optimization of molecular geometry in normal coordinates and use as a tool for simulation of vibrational spectra. The Journal of Chemical Physics, 117(9), 4126. doi:10.1063/1.1498468
-
panther.nmrelaxation.
update_hessian
(grad, grad_old, dx, hessian, update='BFGS')[source]¶ Perform hessian update
Parameters: grad : array_like (N,)
Current gradient
grad_old : array_like (N,)
Previous gradient
dx : array_like (N,)
Step vector x_n - x_(n-1)
hessian : array_like (N, N)
Hessian matrix
update : str
Name of the hessian update to perform, possible values are ‘BFGS’, ‘SR1’ and ‘DFP’
Returns: hessian : array_like
Update hessian matrix
panther.panther module¶
Python package for Anharmonic Thermochemistry
panther.pes module¶
-
panther.pes.
calculate_energies
(images, calc, modes='all')[source]¶ Given a set of images as a nested OrderedDict of Atoms objects and a calculator, calculate the energy for each displaced structure
Parameters: images : OrderedDict
A nested OrderedDict of displaced Atoms objects
calc : calculator instance
ASE calculator
modes : str or list
Mode for which the PES will be calculated
Returns: energies : pandas.DataFrame
DataFrame with the energies per displacement
-
panther.pes.
differentiate
(displacements, energies, order=2)[source]¶ Calculate numerical detivatives using the central difference formula
Parameters: displacements : array_like
energies : DataFrame
order : int
Order of the derivative
Notes
Central difference coefficients taken from [R11]
[R11] Fornberg, B. (1988). Generation of finite difference formulas on arbitrarily spaced grids. Mathematics of Computation, 51(184), 699-699. doi:10.1090/S0025-5718-1988-0935077-0
-
panther.pes.
expandrange
(modestr)[source]¶ Convert a comma separated string of indices and dash separated ranges into a list of integer indices
Parameters: modestr : str Returns: indices : list of ints Examples
>>> from panther.pes import expandrange >>> s = "2,3,5-10,20,25-30" >>> expandrange(s) [2, 3, 5, 6, 7, 8, 9, 10, 20, 25, 26, 27, 28, 29, 30]
-
panther.pes.
fit_potentials
(modeinfo, energies)[source]¶ Fit the potentials with 6th and 4th order polynomials
Parameters: modeinfo : pandas.DataFrame
DataFrame with per mode characteristics, displacements, masses and a flag to mark it a mode is a stretching mode or not
energies : pd.DataFrame
Energies per displacement
Returns: out : (coeffs6o, coeffs4o)
DataFrames with 6th and 4th polynomial coefficients fitted to the potential
panther.plotting module¶
Functions for plotting the each mode and PES fits
panther.vibrations module¶
-
panther.vibrations.
harmonic_vibrational_analysis
(hessian, atoms, proj_translations=True, proj_rotations=False, ascomplex=True, massau=True)[source]¶ Given a force constant matrix (hessian) perform the harmonic vibrational analysis, by calculating the eigevalues and eigenvectors of the mass weighted hessian. Additionally projection of the translational and rotational degrees of freedom can be performed by specifying
proj_translations
andproj_rotations
argsuments.Parameters: hessian : array_like
Force constant (Hessian) matrix in atomic units, should be square and symmetric
atoms : Atoms
ASE atoms object
proj_translations : bool
If
True
translational degrees of freedom will be projected from the hessianproj_rotations : bool
If
True
rotational degrees of freedom will be projected from the hessianmassau : bool
If
True
atomic units of mass will be usedascomplex : bool
If there are complex eigenvalues return the array as complex type otherwise make the complex values negative and return array of reals
Returns: out : (w, v)
Tuple of numpy arrays with hessian square roots of the eigevalues (frequencies) and eiegenvectors in atomic units, both sorted in descending order of eigenvalues
-
panther.vibrations.
project
(atoms, hessian, ndof, proj_translations=True, proj_rotations=False, verbose=False)[source]¶ Project out the translational and/or rotational degrees of freedom from the hessian.
Parameters: atoms : ase.Atoms
Atoms object
ndof : int
Number of degrees of freedom
hessian : array_like
Hessian/force constant matrix
proj_translations : bool
If
True
translational degrees of freedom will be projected from the hessianproj_rotations : bool
If
True
rotational degrees of freedom will be projected from the hessianReturns: proj_hessian : array_like
Hessian matrix with translational and/or rotational degrees of freedom projected out