Molecule

class psi4.core.Molecule

Bases: pybind11_builtins.pybind11_object

Class to store the elements, coordinates, fragmentation pattern, basis sets, charge, multiplicity, etc. of a molecule.

Methods Summary

B787(ref_mol[, do_plot, verbose, atoms_map, …]) Finds shift, rotation, and atom reordering of concern_mol that best aligns with ref_mol.
BFS([seed_atoms, bond_threshold, …]) Detect fragments among real atoms through a breadth-first search (BFS) algorithm.
Z(self, arg0) Nuclear charge of atom
activate_all_fragments(self) Sets all fragments in the molecule to be active
add_atom(self, arg0, arg1, arg2, arg3, arg4, …) Adds to Molecule arg1 an atom with atomic number arg2, Cartesian coordinates in Bohr (arg3, arg4, arg5), atomic symbol arg6, mass arg7, charge arg8 (optional), and lineno arg9 (optional)
atom_at_position(self, arg0, arg1) Tests to see if an atom is at the position arg2 with a given tolerance arg3
basis_on_atom(self, arg0) Gets the label of the orbital basis set on a given atom.
center_of_mass(self) Computes center of mass of molecule (does not translate molecule)
charge(self, arg0) Gets charge of atom
clone(self) Returns a new Molecule identical to arg1
com_fixed(self) Get whether or not COM is fixed
create_molecule_from_string(arg0) Returns a new Molecule with member data from the geometry string arg1 in psi4 format
create_psi4_string_from_molecule(self) Gets a string reexpressing in input format the current states of the molecule
deactivate_all_fragments(self) Sets all fragments in the molecule to be inactive
distance_matrix(self) Returns Matrix of interatom distances
extract_subsets(*args, **kwargs) Overloaded function.
fZ(self, arg0) Nuclear charge of atom arg1 (0-indexed including dummies)
find_point_group(self, arg0) Finds computational molecular point group, user can override this with the symmetry keyword
fix_com(self, arg0) Whether to fix the Cartesian position, or to translate to the C.O.M.
fix_orientation(self, arg0) Fix the orientation at its current frame
form_symmetry_information(self, arg0) Uses the point group object obtain by calling point_group()
from_arrays([geom, elea, elez, elem, mass, …]) Construct Molecule from unvalidated arrays and variables.
from_dict(arg0) Returns a new Molecule constructed from python dictionary.
from_schema(molschema[, return_dict, verbose]) Construct Molecule from non-Psi4 schema.
from_string(molstr[, dtype, name, fix_com, …])
fx(self, arg0) x position of atom arg1 (0-indexed including dummies in Bohr)
fy(self, arg0) y position of atom arg1 (0-indexed including dummies in Bohr)
fz(self, arg0) z position of atom arg1 (0-indexed including dummies in Bohr)
geometry(self) Gets the geometry as a (Natom X 3) matrix of coordinates (in Bohr)
get_fragment_charges(self) Gets the charge of each fragment
get_fragment_multiplicities(self) Gets the multiplicity of each fragment
get_fragment_types(self) Returns a list describing how to handle each fragment {Real, Ghost, Absent}
get_fragments(self) Returns list of pairs of atom ranges defining each fragment from parent molecule(fragments[frag_ind] = <Afirst,Alast+1>)
get_full_point_group(self) Gets point group name such as C3v or S8
get_variable(self, arg0) Checks if variable arg2 is in the list, sets it to val and returns true if it is, and returns false if not
inertia_tensor(self) Returns intertial tensor
input_units_to_au(self) Returns unit conversion to [a0] for geometry
irrep_labels(self) Returns Irreducible Representation symmetry labels
is_variable(self, arg0) Checks if variable arg2 is in the list, returns true if it is, and returns false if not
label(self, arg0) Gets the original label of the atom as given in the input file (C2, H4)
mass(self, atom) Returns mass of atom (0-indexed)
mass_number(self, arg0) Mass number (A) of atom if known, else -1
molecular_charge(self) Gets the molecular charge
move_to_com(self) Moves molecule to center of mass
multiplicity(self) Gets the multiplicity (defined as 2Ms + 1)
nallatom(self) Number of real and dummy atoms
name(self) Gets molecule name
natom(self) Number of real atoms
nfragments(self) Gets the number of fragments in the molecule
nuclear_dipole(*args, **kwargs) Overloaded function.
nuclear_repulsion_energy(self, dipole_field, …) Computes nuclear repulsion energy
nuclear_repulsion_energy_deriv1(self, arg0) Returns first derivative of nuclear repulsion energy as a matrix (natom, 3)
nuclear_repulsion_energy_deriv2(self) Returns second derivative of nuclear repulsion energy as a matrix (natom X 3, natom X 3)
orientation_fixed(self) Get whether or not orientation is fixed
point_group(self) Returns the current point group object
print_bond_angles(self) Print the bond angle geometrical parameters
print_cluster(self) Prints the molecule in Cartesians in input units adding fragment separators
print_distances(self) Print the interatomic distance geometrical parameters
print_in_input_format(self) Prints the molecule as Cartesian or ZMatrix entries, just as inputted.
print_out(self) Prints the molecule in Cartesians in input units to output file
print_out_in_angstrom(self) Prints the molecule in Cartesians in Angstroms to output file
print_out_in_bohr(self) Prints the molecule in Cartesians in Bohr to output file
print_out_of_planes(self) Print the out-of-plane angle geometrical parameters to output file
reinterpret_coordentry(self, arg0) Do reinterpret coordinate entries during update_geometry().
reset_point_group(self, arg0) Overrides symmetry from outside the molecule string
rotational_constants(self) Prints the rotational constants of the molecule
rotational_symmetry_number(self) Returns number of unique orientations of the rigid molecule that only interchange identical atoms
rotor_type(self) Returns rotor type, e.g.
run_dftd3([func, dashlvl, dashparam, …]) Compute dispersion correction using Grimme’s DFTD3 executable.
run_gcp([func, dertype, verbose]) Function to call Grimme’s dftd3 program (http://toc.uni-muenster.de/DFTD3/) to compute the -D correction of level dashlvl using parameters for the functional func.
save_string_xyz(self) Saves the string of an XYZ file to arg2
save_string_xyz_file(self) Saves an XYZ file to arg2
save_xyz_file(self, arg0, arg1) Saves an XYZ file to arg2
schoenflies_symbol(self) Returns the Schoenflies symbol
scramble([do_shift, do_rotate, do_resort, …]) Tester for B787 by shifting, rotating, and atom shuffling ref_mol and checking that the aligner returns the opposite transformation.
set_active_fragment(self, arg0) Sets the specified fragment arg2 to be Real
set_active_fragments(self, arg0) Sets the specified list arg2 of fragments to be Real
set_basis_all_atoms(self, arg0, arg1) Sets basis set arg2 to all atoms
set_basis_by_label(self, arg0, arg1, arg2) Sets basis set arg3 to all atoms with label (e.g., H4) arg2
set_basis_by_symbol(self, arg0, arg1, arg2) Sets basis set arg3 to all atoms with symbol (e.g., H) arg2
set_geometry(self, arg0) Sets the geometry, given a (Natom X 3) matrix arg2 of coordinates (in Bohr)
set_ghost_fragment(self, arg0) Sets the specified fragment arg2 to be Ghost
set_ghost_fragments(self, arg0) Sets the specified list arg2 of fragments to be Ghost
set_input_units_to_au(self, arg0) Sets unit conversion to [a0] for geometry
set_mass(self, atom, mass) Sets mass of atom (0-indexed) to mass
set_molecular_charge(self, arg0) Sets the molecular charge
set_multiplicity(self, arg0) Sets the multiplicity (defined as 2Ms + 1)
set_name(self, arg0) Sets molecule name
set_nuclear_charge(self, arg0, arg1) Set the nuclear charge of the given atom to the value provided.
set_point_group(self, arg0) Sets the molecular point group to the point group object arg2
set_units(self, arg0) Sets units (Angstrom or Bohr) used to define the geometry
set_variable(self, arg0, arg1) Assigns the value arg3 to the variable arg2 in the list of geometry variables, then calls update_geometry()
symbol(self, arg0) Gets the cleaned up label of atom arg2 (C2 => C, H4 = H)
symmetrize(self, arg0) Finds the highest point Abelian point group within the specified tolerance, and forces the geometry to have that symmetry.
symmetry_from_input(self) Returns the symmetry specified in the input
to_arrays() Exports coordinate info into NumPy arrays.
to_dict([force_c1, force_units, np_out]) Serializes instance into Molecule dictionary.
to_schema(dtype[, units, return_type]) Serializes instance into JSON or YAML according to schema dtype.
to_string(dtype[, units, atom_format, …]) Format a string representation of QM molecule.
translate(self, arg0) Translates molecule by arg2
units(self) Returns units used to define the geometry, i.e.
update_geometry(self) Reevaluates the geometry with current variable values, orientation directives, etc.
x(self, arg0) x position of atom
y(self, arg0) y position of atom
z(self, arg0) z position of atom

Methods Documentation

B787(ref_mol, do_plot=False, verbose=1, atoms_map=False, run_resorting=False, mols_align=False, run_to_completion=False, uno_cutoff=0.001, run_mirror=False)

Finds shift, rotation, and atom reordering of concern_mol that best aligns with ref_mol.

Wraps qcdb.align.B787() for qcdb.Molecule or psi4.core.Molecule. Employs the Kabsch, Hungarian, and Uno algorithms to exhaustively locate the best alignment for non-oriented, non-ordered structures.

Parameters:
  • concern_mol (qcdb.Molecule or psi4.core.Molecule) – Molecule of concern, to be shifted, rotated, and reordered into best coincidence with ref_mol.
  • ref_mol (qcdb.Molecule or psi4.core.Molecule) – Molecule to match.
  • atoms_map (bool, optional) – Whether atom1 of ref_mol corresponds to atom1 of concern_mol, etc. If true, specifying True can save much time.
  • mols_align (bool, optional) – Whether ref_mol and concern_mol have identical geometries by eye (barring orientation or atom mapping) and expected final RMSD = 0. If True, procedure is truncated when RMSD condition met, saving time.
  • do_plot (bool, optional) – Pops up a mpl plot showing before, after, and ref geometries.
  • run_to_completion (bool, optional) – Run reorderings to completion (past RMSD = 0) even if unnecessary because mols_align=True. Used to test worst-case timings.
  • run_resorting (bool, optional) – Run the resorting machinery even if unnecessary because atoms_map=True.
  • uno_cutoff (float, optional) – TODO
  • run_mirror (bool, optional) – Run alternate geometries potentially allowing best match to ref_mol from mirror image of concern_mol. Only run if system confirmed to be nonsuperimposable upon mirror reflection.
Returns:

First item is RMSD [A] between ref_mol and the optimally aligned geometry computed. Second item is a AlignmentMill namedtuple with fields (shift, rotation, atommap, mirror) that prescribe the transformation from concern_mol and the optimally aligned geometry. Third item is a crude charge-, multiplicity-, fragment-less Molecule at optimally aligned (and atom-ordered) geometry. Return type determined by concern_mol type.

Return type:

float, tuple, qcdb.Molecule or psi4.core.Molecule

BFS(seed_atoms=None, bond_threshold=1.2, return_arrays=False, return_molecules=False, return_molecule=False)

Detect fragments among real atoms through a breadth-first search (BFS) algorithm.

Parameters:
  • self (qcdb.Molecule or psi4.core.Molecule) –
  • seed_atoms (list, optional) – List of lists of atoms (0-indexed) belonging to independent fragments. Useful to prompt algorithm or to define intramolecular fragments through border atoms. Example: [[1, 0], [2]]
  • bond_threshold (float, optional) – Factor beyond average of covalent radii to determine bond cutoff.
  • return_arrays (bool, optional) – If True, also return fragments as list of arrays.
  • return_molecules (bool, optional) – If True, also return fragments as list of Molecules.
  • return_molecule (bool, optional) – If True, also return one big Molecule with fragmentation encoded.
Returns:

  • bfs_map (list of lists) – Array of atom indices (0-indexed) of detected fragments.
  • bfs_arrays (tuple of lists of ndarray, optional) – geom, mass, elem info per-fragment. Only provided if return_arrays is True.
  • bfs_molecules (list of qcdb.Molecule or psi4.core.Molecule, optional) – List of molecules, each built from one fragment. Center and orientation of fragments is fixed so orientation info from self is not lost. Loses chgmult and ghost/dummy info from self and contains default chgmult. Only provided if return_molecules is True. Returned are of same type as self.
  • bfs_molecule (qcdb.Molecule or psi4.core.Molecule, optional) – Single molecule with same number of real atoms as self with atoms reordered into adjacent fragments and fragment markers inserted. Loses ghost/dummy info from self; keeps total charge but not total mult. Only provided if return_molecule is True. Returned is of same type as self.

Notes

Relies upon van der Waals radii and so faulty for close (especially
hydrogen-bonded) fragments. See seed_atoms.

Any existing fragmentation info/chgmult encoded in self is lost.

Original code from Michael S. Marshall, linear-scaling algorithm from Trent M. Parker, revamped by Lori A. Burns

Z(self: psi4.core.Molecule, arg0: int) → float

Nuclear charge of atom

activate_all_fragments(self: psi4.core.Molecule) → None

Sets all fragments in the molecule to be active

add_atom(self: psi4.core.Molecule, arg0: float, arg1: float, arg2: float, arg3: float, arg4: str, arg5: float, arg6: float, arg7: str, arg8: int) → None

Adds to Molecule arg1 an atom with atomic number arg2, Cartesian coordinates in Bohr (arg3, arg4, arg5), atomic symbol arg6, mass arg7, charge arg8 (optional), and lineno arg9 (optional)

atom_at_position(self: psi4.core.Molecule, arg0: float, arg1: float) → int

Tests to see if an atom is at the position arg2 with a given tolerance arg3

basis_on_atom(self: psi4.core.Molecule, arg0: int) → str

Gets the label of the orbital basis set on a given atom.

center_of_mass(self: psi4.core.Molecule) → psi4.core.Vector3

Computes center of mass of molecule (does not translate molecule)

charge(self: psi4.core.Molecule, arg0: int) → float

Gets charge of atom

clone(self: psi4.core.Molecule) → psi4.core.Molecule

Returns a new Molecule identical to arg1

com_fixed(self: psi4.core.Molecule) → bool

Get whether or not COM is fixed

create_molecule_from_string(arg0: str) → psi4.core.Molecule

Returns a new Molecule with member data from the geometry string arg1 in psi4 format

create_psi4_string_from_molecule(self: psi4.core.Molecule) → str

Gets a string reexpressing in input format the current states of the molecule

deactivate_all_fragments(self: psi4.core.Molecule) → None

Sets all fragments in the molecule to be inactive

distance_matrix(self: psi4.core.Molecule) → psi4.core.Matrix

Returns Matrix of interatom distances

extract_subsets(*args, **kwargs)

Overloaded function.

  1. extract_subsets(self: psi4.core.Molecule, arg0: List[int], arg1: List[int]) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragments Real and arg3 fragments Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: List[int], arg1: int) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragments Real and arg3 fragment Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: int, arg1: List[int]) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragment Real and arg3 fragments Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: int, arg1: int) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragment Real and arg3 fragment Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: List[int]) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragments Real

  1. extract_subsets(self: psi4.core.Molecule, arg0: int) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragment Real

fZ(self: psi4.core.Molecule, arg0: int) → float

Nuclear charge of atom arg1 (0-indexed including dummies)

find_point_group(self: psi4.core.Molecule, arg0: float) → psi4.core.PointGroup

Finds computational molecular point group, user can override this with the symmetry keyword

fix_com(self: psi4.core.Molecule, arg0: bool) → None

Whether to fix the Cartesian position, or to translate to the C.O.M.

fix_orientation(self: psi4.core.Molecule, arg0: bool) → None

Fix the orientation at its current frame

form_symmetry_information(self: psi4.core.Molecule, arg0: float) → None

Uses the point group object obtain by calling point_group()

classmethod from_arrays(geom=None, elea=None, elez=None, elem=None, mass=None, real=None, elbl=None, name=None, units='Angstrom', input_units_to_au=None, fix_com=None, fix_orientation=None, fix_symmetry=None, fragment_separators=None, fragment_charges=None, fragment_multiplicities=None, molecular_charge=None, molecular_multiplicity=None, missing_enabled_return='error', tooclose=0.1, zero_ghost_fragments=False, nonphysical=False, mtol=0.001, verbose=1, return_dict=False)

Construct Molecule from unvalidated arrays and variables.

Light wrapper around from_arrays() that is a full-featured constructor to dictionary representa- tion of Molecule. This follows one step further to return Molecule instance.

:param See from_arrays().:

Returns:
Return type:psi4.core.Molecule
from_dict(arg0: dict) → psi4.core.Molecule

Returns a new Molecule constructed from python dictionary. In progress: name and capabilities should not be relied upon

classmethod from_schema(molschema, return_dict=False, verbose=1)

Construct Molecule from non-Psi4 schema.

Light wrapper around from_arrays().

Parameters:
  • molschema (dict) – Dictionary form of Molecule following known schema.
  • return_dict (bool, optional) – Additionally return Molecule dictionary intermediate.
  • verbose (int, optional) – Amount of printing.
Returns:

  • mol (psi4.core.Molecule)
  • molrec (dict, optional) – Dictionary representation of instance. Only provided if return_dict is True.

classmethod from_string(molstr, dtype=None, name=None, fix_com=None, fix_orientation=None, fix_symmetry=None, return_dict=False, enable_qm=True, enable_efp=True, missing_enabled_return_qm='none', missing_enabled_return_efp='none', verbose=1)
fx(self: psi4.core.Molecule, arg0: int) → float

x position of atom arg1 (0-indexed including dummies in Bohr)

fy(self: psi4.core.Molecule, arg0: int) → float

y position of atom arg1 (0-indexed including dummies in Bohr)

fz(self: psi4.core.Molecule, arg0: int) → float

z position of atom arg1 (0-indexed including dummies in Bohr)

geometry(self: psi4.core.Molecule) → psi4.core.Matrix

Gets the geometry as a (Natom X 3) matrix of coordinates (in Bohr)

get_fragment_charges(self: psi4.core.Molecule) → List[int]

Gets the charge of each fragment

get_fragment_multiplicities(self: psi4.core.Molecule) → List[int]

Gets the multiplicity of each fragment

get_fragment_types(self: psi4.core.Molecule) → List[str]

Returns a list describing how to handle each fragment {Real, Ghost, Absent}

get_fragments(self: psi4.core.Molecule) → List[Tuple[int, int]]

Returns list of pairs of atom ranges defining each fragment from parent molecule(fragments[frag_ind] = <Afirst,Alast+1>)

get_full_point_group(self: psi4.core.Molecule) → str

Gets point group name such as C3v or S8

get_variable(self: psi4.core.Molecule, arg0: str) → float

Checks if variable arg2 is in the list, sets it to val and returns true if it is, and returns false if not

inertia_tensor(self: psi4.core.Molecule) → psi4.core.Matrix

Returns intertial tensor

input_units_to_au(self: psi4.core.Molecule) → float

Returns unit conversion to [a0] for geometry

irrep_labels(self: psi4.core.Molecule) → List[str]

Returns Irreducible Representation symmetry labels

is_variable(self: psi4.core.Molecule, arg0: str) → bool

Checks if variable arg2 is in the list, returns true if it is, and returns false if not

label(self: psi4.core.Molecule, arg0: int) → str

Gets the original label of the atom as given in the input file (C2, H4)

mass(self: psi4.core.Molecule, atom: int) → float

Returns mass of atom (0-indexed)

mass_number(self: psi4.core.Molecule, arg0: int) → int

Mass number (A) of atom if known, else -1

molecular_charge(self: psi4.core.Molecule) → int

Gets the molecular charge

move_to_com(self: psi4.core.Molecule) → None

Moves molecule to center of mass

multiplicity(self: psi4.core.Molecule) → int

Gets the multiplicity (defined as 2Ms + 1)

nallatom(self: psi4.core.Molecule) → int

Number of real and dummy atoms

name(self: psi4.core.Molecule) → str

Gets molecule name

natom(self: psi4.core.Molecule) → int

Number of real atoms

nfragments(self: psi4.core.Molecule) → int

Gets the number of fragments in the molecule

nuclear_dipole(*args, **kwargs)

Overloaded function.

  1. nuclear_dipole(self: psi4.core.Molecule, arg0: psi4.core.Vector3) -> psi4.core.Vector3

Gets the nuclear contribution to the dipole, withe respect to a specified origin

  1. nuclear_dipole(self: psi4.core.Molecule) -> psi4.core.Vector3

Gets the nuclear contribution to the dipole, withe respect to the origin

nuclear_repulsion_energy(self: psi4.core.Molecule, dipole_field: List[float[3]]=[0.0, 0.0, 0.0]) → float

Computes nuclear repulsion energy

nuclear_repulsion_energy_deriv1(self: psi4.core.Molecule, arg0: List[float[3]]) → psi4.core.Matrix

Returns first derivative of nuclear repulsion energy as a matrix (natom, 3)

nuclear_repulsion_energy_deriv2(self: psi4.core.Molecule) → psi4.core.Matrix

Returns second derivative of nuclear repulsion energy as a matrix (natom X 3, natom X 3)

orientation_fixed(self: psi4.core.Molecule) → bool

Get whether or not orientation is fixed

point_group(self: psi4.core.Molecule) → psi4.core.PointGroup

Returns the current point group object

print_bond_angles(self: psi4.core.Molecule) → None

Print the bond angle geometrical parameters

print_cluster(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in input units adding fragment separators

print_distances(self: psi4.core.Molecule) → None

Print the interatomic distance geometrical parameters

print_in_input_format(self: psi4.core.Molecule) → None

Prints the molecule as Cartesian or ZMatrix entries, just as inputted.

print_out(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in input units to output file

print_out_in_angstrom(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in Angstroms to output file

print_out_in_bohr(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in Bohr to output file

print_out_of_planes(self: psi4.core.Molecule) → None

Print the out-of-plane angle geometrical parameters to output file

reinterpret_coordentry(self: psi4.core.Molecule, arg0: bool) → None

Do reinterpret coordinate entries during update_geometry().

reset_point_group(self: psi4.core.Molecule, arg0: str) → None

Overrides symmetry from outside the molecule string

rotational_constants(self: psi4.core.Molecule) → psi4.core.Vector

Prints the rotational constants of the molecule

rotational_symmetry_number(self: psi4.core.Molecule) → int

Returns number of unique orientations of the rigid molecule that only interchange identical atoms

rotor_type(self: psi4.core.Molecule) → str

Returns rotor type, e.g. ‘RT_ATOM’ or ‘RT_SYMMETRIC_TOP’

run_dftd3(func=None, dashlvl=None, dashparam=None, dertype=None, verbose=False)

Compute dispersion correction using Grimme’s DFTD3 executable.

Function to call Grimme’s dftd3 program to compute the -D correction of level dashlvl using parameters for the functional func. dashparam can supply a full set of dispersion parameters in the absence of func or individual overrides in the presence of func.

The DFTD3 executable must be independently compiled and found in PATH or PSIPATH.

Parameters:
  • mol (qcdb.Molecule or psi4.core.Molecule or str) – Molecule on which to run dispersion calculation. Both qcdb and psi4.core Molecule classes have been extended by this method, so either allowed. Alternately, a string that can be instantiated into a qcdb.Molecule.
  • func (str or None) – Density functional (Psi4, not Turbomole, names) for which to load parameters from dashcoeff[dashlvl][func]. This is not passed to DFTD3 and thus may be a dummy or None. Any or all parameters initialized can be overwritten via dashparam.
  • dashlvl ({'d2p4', 'd2gr', 'd3zero', 'd3bj', 'd3mzero', d3mbj', 'd', 'd2', 'd3', 'd3m'}) – Flavor of a posteriori dispersion correction for which to load parameters and call procedure in DFTD3. Must be a keys in dashcoeff dict (or a key in dashalias that resolves to one).
  • dashparam (dict, optional) – Dictionary of the same keys as dashcoeff[dashlvl] used to override any or all values initialized by dashcoeff[dashlvl][func].
  • dertype ({None, 0, 'none', 'energy', 1, 'first', 'gradient'}, optional) – Maximum derivative level at which to run DFTD3. For large mol, energy-only calculations can be significantly more efficient. Also controls return values, see below.
  • verbose (bool, optional) – When True, additionally include DFTD3 output in output.
Returns:

  • energy (float, optional) – When dertype is 0, energy [Eh].
  • gradient (list of lists of floats or psi4.core.Matrix, optional) – When dertype is 1, (nat, 3) gradient [Eh/a0].
  • (energy, gradient) (float and list of lists of floats or psi4.core.Matrix, optional) – When dertype is unspecified, both energy [Eh] and (nat, 3) gradient [Eh/a0].

Notes

research site: https://www.chemie.uni-bonn.de/pctc/mulliken-center/software/dft-d3 Psi4 mode: When psi4 the python module is importable at import qcdb

time, Psi4 mode is activated, with the following alterations: * output goes to output file * gradient returned as psi4.core.Matrix, not list o’lists * scratch is written to randomly named subdirectory of psi scratch * psivar “DISPERSION CORRECTION ENERGY” is set * verbose triggered when PRINT keywork of SCF module >=3
run_gcp(func=None, dertype=None, verbose=False)

Function to call Grimme’s dftd3 program (http://toc.uni-muenster.de/DFTD3/) to compute the -D correction of level dashlvl using parameters for the functional func. The dictionary dashparam can be used to supply a full set of dispersion parameters in the absense of func or to supply individual overrides in the presence of func. Returns energy if dertype is 0, gradient if dertype is 1, else tuple of energy and gradient if dertype unspecified. The dftd3 executable must be independently compiled and found in PATH or PSIPATH. self may be either a qcdb.Molecule (sensibly) or a psi4.Molecule (works b/c psi4.Molecule has been extended by this method py-side and only public interface fns used) or a string that can be instantiated into a qcdb.Molecule.

save_string_xyz(self: psi4.core.Molecule) → str

Saves the string of an XYZ file to arg2

save_string_xyz_file(self: psi4.core.Molecule) → str

Saves an XYZ file to arg2

save_xyz_file(self: psi4.core.Molecule, arg0: str, arg1: bool) → None

Saves an XYZ file to arg2

schoenflies_symbol(self: psi4.core.Molecule) → str

Returns the Schoenflies symbol

scramble(do_shift=True, do_rotate=True, do_resort=True, deflection=1.0, do_mirror=False, do_plot=False, run_to_completion=False, run_resorting=False, verbose=1)

Tester for B787 by shifting, rotating, and atom shuffling ref_mol and checking that the aligner returns the opposite transformation.

Parameters:
  • ref_mol (qcdb.Molecule or psi4.core.Molecule) – Molecule to perturb.
  • do_shift (bool or array-like, optional) – Whether to generate a random atom shift on interval [-3, 3) in each dimension (True) or leave at current origin. To shift by a specified vector, supply a 3-element list.
  • do_rotate (bool or array-like, optional) – Whether to generate a random 3D rotation according to algorithm of Arvo. To rotate by a specified matrix, supply a 9-element list of lists.
  • do_resort (bool or array-like, optional) – Whether to shuffle atoms (True) or leave 1st atom 1st, etc. (False). To specify shuffle, supply a nat-element list of indices.
  • deflection (float, optional) – If do_rotate, how random a rotation: 0.0 is no change, 0.1 is small perturbation, 1.0 is completely random.
  • do_mirror (bool, optional) – Whether to construct the mirror image structure by inverting y-axis.
  • do_plot (bool, optional) – Pops up a mpl plot showing before, after, and ref geometries.
  • run_to_completion (bool, optional) – By construction, scrambled systems are fully alignable (final RMSD=0). Even so, True turns off the mechanism to stop when RMSD reaches zero and instead proceed to worst possible time.
  • run_resorting (bool, optional) – Even if atoms not shuffled, test the resorting machinery.
  • verbose (int, optional) – Print level.
Returns:

Return type:

None

set_active_fragment(self: psi4.core.Molecule, arg0: int) → None

Sets the specified fragment arg2 to be Real

set_active_fragments(self: psi4.core.Molecule, arg0: List[int]) → None

Sets the specified list arg2 of fragments to be Real

set_basis_all_atoms(self: psi4.core.Molecule, arg0: str, arg1: str) → None

Sets basis set arg2 to all atoms

set_basis_by_label(self: psi4.core.Molecule, arg0: str, arg1: str, arg2: str) → None

Sets basis set arg3 to all atoms with label (e.g., H4) arg2

set_basis_by_symbol(self: psi4.core.Molecule, arg0: str, arg1: str, arg2: str) → None

Sets basis set arg3 to all atoms with symbol (e.g., H) arg2

set_geometry(self: psi4.core.Molecule, arg0: psi4.core.Matrix) → None

Sets the geometry, given a (Natom X 3) matrix arg2 of coordinates (in Bohr)

set_ghost_fragment(self: psi4.core.Molecule, arg0: int) → None

Sets the specified fragment arg2 to be Ghost

set_ghost_fragments(self: psi4.core.Molecule, arg0: List[int]) → None

Sets the specified list arg2 of fragments to be Ghost

set_input_units_to_au(self: psi4.core.Molecule, arg0: float) → None

Sets unit conversion to [a0] for geometry

set_mass(self: psi4.core.Molecule, atom: int, mass: float) → None

Sets mass of atom (0-indexed) to mass

set_molecular_charge(self: psi4.core.Molecule, arg0: int) → None

Sets the molecular charge

set_multiplicity(self: psi4.core.Molecule, arg0: int) → None

Sets the multiplicity (defined as 2Ms + 1)

set_name(self: psi4.core.Molecule, arg0: str) → None

Sets molecule name

set_nuclear_charge(self: psi4.core.Molecule, arg0: int, arg1: float) → None

Set the nuclear charge of the given atom to the value provided.

set_point_group(self: psi4.core.Molecule, arg0: psi4.core.PointGroup) → None

Sets the molecular point group to the point group object arg2

set_units(self: psi4.core.Molecule, arg0: psi4.core.GeometryUnits) → None

Sets units (Angstrom or Bohr) used to define the geometry

set_variable(self: psi4.core.Molecule, arg0: str, arg1: float) → None

Assigns the value arg3 to the variable arg2 in the list of geometry variables, then calls update_geometry()

symbol(self: psi4.core.Molecule, arg0: int) → str

Gets the cleaned up label of atom arg2 (C2 => C, H4 = H)

symmetrize(self: psi4.core.Molecule, arg0: float) → None

Finds the highest point Abelian point group within the specified tolerance, and forces the geometry to have that symmetry.

symmetry_from_input(self: psi4.core.Molecule) → str

Returns the symmetry specified in the input

to_arrays()

Exports coordinate info into NumPy arrays.

Returns:
  • geom, mass, elem, elez, uniq (ndarray, ndarray, ndarray, ndarray, ndarray) – (nat, 3) geometry [a0]. (nat,) mass [u]. (nat,) element symbol. (nat,) atomic number. (nat,) hash of element symbol and mass. Note that coordinate, orientation, and element information is preserved but fragmentation, chgmult, and dummy/ghost is lost.
  • Usage
  • —–
  • geom, mass, elem, elez, uniq = molinstance.to_arrays()
to_dict(force_c1=False, force_units=False, np_out=True)

Serializes instance into Molecule dictionary.

to_schema(dtype, units='Angstrom', return_type='json')

Serializes instance into JSON or YAML according to schema dtype.

to_string(dtype, units='Angstrom', atom_format=None, ghost_format=None, width=17, prec=12)

Format a string representation of QM molecule.

translate(self: psi4.core.Molecule, arg0: psi4.core.Vector3) → None

Translates molecule by arg2

units(self: psi4.core.Molecule) → str

Returns units used to define the geometry, i.e. ‘Angstrom’ or ‘Bohr’

update_geometry(self: psi4.core.Molecule) → None

Reevaluates the geometry with current variable values, orientation directives, etc. Must be called after initial Molecule definition by string.

x(self: psi4.core.Molecule, arg0: int) → float

x position of atom

y(self: psi4.core.Molecule, arg0: int) → float

y position of atom

z(self: psi4.core.Molecule, arg0: int) → float

z position of atom

B787(ref_mol, do_plot=False, verbose=1, atoms_map=False, run_resorting=False, mols_align=False, run_to_completion=False, uno_cutoff=0.001, run_mirror=False)

Finds shift, rotation, and atom reordering of concern_mol that best aligns with ref_mol.

Wraps qcdb.align.B787() for qcdb.Molecule or psi4.core.Molecule. Employs the Kabsch, Hungarian, and Uno algorithms to exhaustively locate the best alignment for non-oriented, non-ordered structures.

Parameters:
  • concern_mol (qcdb.Molecule or psi4.core.Molecule) – Molecule of concern, to be shifted, rotated, and reordered into best coincidence with ref_mol.
  • ref_mol (qcdb.Molecule or psi4.core.Molecule) – Molecule to match.
  • atoms_map (bool, optional) – Whether atom1 of ref_mol corresponds to atom1 of concern_mol, etc. If true, specifying True can save much time.
  • mols_align (bool, optional) – Whether ref_mol and concern_mol have identical geometries by eye (barring orientation or atom mapping) and expected final RMSD = 0. If True, procedure is truncated when RMSD condition met, saving time.
  • do_plot (bool, optional) – Pops up a mpl plot showing before, after, and ref geometries.
  • run_to_completion (bool, optional) – Run reorderings to completion (past RMSD = 0) even if unnecessary because mols_align=True. Used to test worst-case timings.
  • run_resorting (bool, optional) – Run the resorting machinery even if unnecessary because atoms_map=True.
  • uno_cutoff (float, optional) – TODO
  • run_mirror (bool, optional) – Run alternate geometries potentially allowing best match to ref_mol from mirror image of concern_mol. Only run if system confirmed to be nonsuperimposable upon mirror reflection.
Returns:

First item is RMSD [A] between ref_mol and the optimally aligned geometry computed. Second item is a AlignmentMill namedtuple with fields (shift, rotation, atommap, mirror) that prescribe the transformation from concern_mol and the optimally aligned geometry. Third item is a crude charge-, multiplicity-, fragment-less Molecule at optimally aligned (and atom-ordered) geometry. Return type determined by concern_mol type.

Return type:

float, tuple, qcdb.Molecule or psi4.core.Molecule

BFS(seed_atoms=None, bond_threshold=1.2, return_arrays=False, return_molecules=False, return_molecule=False)

Detect fragments among real atoms through a breadth-first search (BFS) algorithm.

Parameters:
  • self (qcdb.Molecule or psi4.core.Molecule) –
  • seed_atoms (list, optional) – List of lists of atoms (0-indexed) belonging to independent fragments. Useful to prompt algorithm or to define intramolecular fragments through border atoms. Example: [[1, 0], [2]]
  • bond_threshold (float, optional) – Factor beyond average of covalent radii to determine bond cutoff.
  • return_arrays (bool, optional) – If True, also return fragments as list of arrays.
  • return_molecules (bool, optional) – If True, also return fragments as list of Molecules.
  • return_molecule (bool, optional) – If True, also return one big Molecule with fragmentation encoded.
Returns:

  • bfs_map (list of lists) – Array of atom indices (0-indexed) of detected fragments.
  • bfs_arrays (tuple of lists of ndarray, optional) – geom, mass, elem info per-fragment. Only provided if return_arrays is True.
  • bfs_molecules (list of qcdb.Molecule or psi4.core.Molecule, optional) – List of molecules, each built from one fragment. Center and orientation of fragments is fixed so orientation info from self is not lost. Loses chgmult and ghost/dummy info from self and contains default chgmult. Only provided if return_molecules is True. Returned are of same type as self.
  • bfs_molecule (qcdb.Molecule or psi4.core.Molecule, optional) – Single molecule with same number of real atoms as self with atoms reordered into adjacent fragments and fragment markers inserted. Loses ghost/dummy info from self; keeps total charge but not total mult. Only provided if return_molecule is True. Returned is of same type as self.

Notes

Relies upon van der Waals radii and so faulty for close (especially
hydrogen-bonded) fragments. See seed_atoms.

Any existing fragmentation info/chgmult encoded in self is lost.

Original code from Michael S. Marshall, linear-scaling algorithm from Trent M. Parker, revamped by Lori A. Burns

Z(self: psi4.core.Molecule, arg0: int) → float

Nuclear charge of atom

activate_all_fragments(self: psi4.core.Molecule) → None

Sets all fragments in the molecule to be active

add_atom(self: psi4.core.Molecule, arg0: float, arg1: float, arg2: float, arg3: float, arg4: str, arg5: float, arg6: float, arg7: str, arg8: int) → None

Adds to Molecule arg1 an atom with atomic number arg2, Cartesian coordinates in Bohr (arg3, arg4, arg5), atomic symbol arg6, mass arg7, charge arg8 (optional), and lineno arg9 (optional)

atom_at_position(self: psi4.core.Molecule, arg0: float, arg1: float) → int

Tests to see if an atom is at the position arg2 with a given tolerance arg3

basis_on_atom(self: psi4.core.Molecule, arg0: int) → str

Gets the label of the orbital basis set on a given atom.

center_of_mass(self: psi4.core.Molecule) → psi4.core.Vector3

Computes center of mass of molecule (does not translate molecule)

charge(self: psi4.core.Molecule, arg0: int) → float

Gets charge of atom

clone(self: psi4.core.Molecule) → psi4.core.Molecule

Returns a new Molecule identical to arg1

com_fixed(self: psi4.core.Molecule) → bool

Get whether or not COM is fixed

create_molecule_from_string(arg0: str) → psi4.core.Molecule

Returns a new Molecule with member data from the geometry string arg1 in psi4 format

create_psi4_string_from_molecule(self: psi4.core.Molecule) → str

Gets a string reexpressing in input format the current states of the molecule

deactivate_all_fragments(self: psi4.core.Molecule) → None

Sets all fragments in the molecule to be inactive

distance_matrix(self: psi4.core.Molecule) → psi4.core.Matrix

Returns Matrix of interatom distances

extract_subsets(*args, **kwargs)

Overloaded function.

  1. extract_subsets(self: psi4.core.Molecule, arg0: List[int], arg1: List[int]) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragments Real and arg3 fragments Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: List[int], arg1: int) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragments Real and arg3 fragment Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: int, arg1: List[int]) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragment Real and arg3 fragments Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: int, arg1: int) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragment Real and arg3 fragment Ghost

  1. extract_subsets(self: psi4.core.Molecule, arg0: List[int]) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragments Real

  1. extract_subsets(self: psi4.core.Molecule, arg0: int) -> psi4.core.Molecule

Returns copy of arg1 with arg2 fragment Real

fZ(self: psi4.core.Molecule, arg0: int) → float

Nuclear charge of atom arg1 (0-indexed including dummies)

find_point_group(self: psi4.core.Molecule, arg0: float) → psi4.core.PointGroup

Finds computational molecular point group, user can override this with the symmetry keyword

fix_com(self: psi4.core.Molecule, arg0: bool) → None

Whether to fix the Cartesian position, or to translate to the C.O.M.

fix_orientation(self: psi4.core.Molecule, arg0: bool) → None

Fix the orientation at its current frame

form_symmetry_information(self: psi4.core.Molecule, arg0: float) → None

Uses the point group object obtain by calling point_group()

classmethod from_arrays(geom=None, elea=None, elez=None, elem=None, mass=None, real=None, elbl=None, name=None, units='Angstrom', input_units_to_au=None, fix_com=None, fix_orientation=None, fix_symmetry=None, fragment_separators=None, fragment_charges=None, fragment_multiplicities=None, molecular_charge=None, molecular_multiplicity=None, missing_enabled_return='error', tooclose=0.1, zero_ghost_fragments=False, nonphysical=False, mtol=0.001, verbose=1, return_dict=False)

Construct Molecule from unvalidated arrays and variables.

Light wrapper around from_arrays() that is a full-featured constructor to dictionary representa- tion of Molecule. This follows one step further to return Molecule instance.

:param See from_arrays().:

Returns:
Return type:psi4.core.Molecule
from_dict(arg0: dict) → psi4.core.Molecule

Returns a new Molecule constructed from python dictionary. In progress: name and capabilities should not be relied upon

classmethod from_schema(molschema, return_dict=False, verbose=1)

Construct Molecule from non-Psi4 schema.

Light wrapper around from_arrays().

Parameters:
  • molschema (dict) – Dictionary form of Molecule following known schema.
  • return_dict (bool, optional) – Additionally return Molecule dictionary intermediate.
  • verbose (int, optional) – Amount of printing.
Returns:

  • mol (psi4.core.Molecule)
  • molrec (dict, optional) – Dictionary representation of instance. Only provided if return_dict is True.

classmethod from_string(molstr, dtype=None, name=None, fix_com=None, fix_orientation=None, fix_symmetry=None, return_dict=False, enable_qm=True, enable_efp=True, missing_enabled_return_qm='none', missing_enabled_return_efp='none', verbose=1)
fx(self: psi4.core.Molecule, arg0: int) → float

x position of atom arg1 (0-indexed including dummies in Bohr)

fy(self: psi4.core.Molecule, arg0: int) → float

y position of atom arg1 (0-indexed including dummies in Bohr)

fz(self: psi4.core.Molecule, arg0: int) → float

z position of atom arg1 (0-indexed including dummies in Bohr)

geometry(self: psi4.core.Molecule) → psi4.core.Matrix

Gets the geometry as a (Natom X 3) matrix of coordinates (in Bohr)

get_fragment_charges(self: psi4.core.Molecule) → List[int]

Gets the charge of each fragment

get_fragment_multiplicities(self: psi4.core.Molecule) → List[int]

Gets the multiplicity of each fragment

get_fragment_types(self: psi4.core.Molecule) → List[str]

Returns a list describing how to handle each fragment {Real, Ghost, Absent}

get_fragments(self: psi4.core.Molecule) → List[Tuple[int, int]]

Returns list of pairs of atom ranges defining each fragment from parent molecule(fragments[frag_ind] = <Afirst,Alast+1>)

get_full_point_group(self: psi4.core.Molecule) → str

Gets point group name such as C3v or S8

get_variable(self: psi4.core.Molecule, arg0: str) → float

Checks if variable arg2 is in the list, sets it to val and returns true if it is, and returns false if not

inertia_tensor(self: psi4.core.Molecule) → psi4.core.Matrix

Returns intertial tensor

input_units_to_au(self: psi4.core.Molecule) → float

Returns unit conversion to [a0] for geometry

irrep_labels(self: psi4.core.Molecule) → List[str]

Returns Irreducible Representation symmetry labels

is_variable(self: psi4.core.Molecule, arg0: str) → bool

Checks if variable arg2 is in the list, returns true if it is, and returns false if not

label(self: psi4.core.Molecule, arg0: int) → str

Gets the original label of the atom as given in the input file (C2, H4)

mass(self: psi4.core.Molecule, atom: int) → float

Returns mass of atom (0-indexed)

mass_number(self: psi4.core.Molecule, arg0: int) → int

Mass number (A) of atom if known, else -1

molecular_charge(self: psi4.core.Molecule) → int

Gets the molecular charge

move_to_com(self: psi4.core.Molecule) → None

Moves molecule to center of mass

multiplicity(self: psi4.core.Molecule) → int

Gets the multiplicity (defined as 2Ms + 1)

nallatom(self: psi4.core.Molecule) → int

Number of real and dummy atoms

name(self: psi4.core.Molecule) → str

Gets molecule name

natom(self: psi4.core.Molecule) → int

Number of real atoms

nfragments(self: psi4.core.Molecule) → int

Gets the number of fragments in the molecule

nuclear_dipole(*args, **kwargs)

Overloaded function.

  1. nuclear_dipole(self: psi4.core.Molecule, arg0: psi4.core.Vector3) -> psi4.core.Vector3

Gets the nuclear contribution to the dipole, withe respect to a specified origin

  1. nuclear_dipole(self: psi4.core.Molecule) -> psi4.core.Vector3

Gets the nuclear contribution to the dipole, withe respect to the origin

nuclear_repulsion_energy(self: psi4.core.Molecule, dipole_field: List[float[3]]=[0.0, 0.0, 0.0]) → float

Computes nuclear repulsion energy

nuclear_repulsion_energy_deriv1(self: psi4.core.Molecule, arg0: List[float[3]]) → psi4.core.Matrix

Returns first derivative of nuclear repulsion energy as a matrix (natom, 3)

nuclear_repulsion_energy_deriv2(self: psi4.core.Molecule) → psi4.core.Matrix

Returns second derivative of nuclear repulsion energy as a matrix (natom X 3, natom X 3)

orientation_fixed(self: psi4.core.Molecule) → bool

Get whether or not orientation is fixed

point_group(self: psi4.core.Molecule) → psi4.core.PointGroup

Returns the current point group object

print_bond_angles(self: psi4.core.Molecule) → None

Print the bond angle geometrical parameters

print_cluster(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in input units adding fragment separators

print_distances(self: psi4.core.Molecule) → None

Print the interatomic distance geometrical parameters

print_in_input_format(self: psi4.core.Molecule) → None

Prints the molecule as Cartesian or ZMatrix entries, just as inputted.

print_out(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in input units to output file

print_out_in_angstrom(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in Angstroms to output file

print_out_in_bohr(self: psi4.core.Molecule) → None

Prints the molecule in Cartesians in Bohr to output file

print_out_of_planes(self: psi4.core.Molecule) → None

Print the out-of-plane angle geometrical parameters to output file

reinterpret_coordentry(self: psi4.core.Molecule, arg0: bool) → None

Do reinterpret coordinate entries during update_geometry().

reset_point_group(self: psi4.core.Molecule, arg0: str) → None

Overrides symmetry from outside the molecule string

rotational_constants(self: psi4.core.Molecule) → psi4.core.Vector

Prints the rotational constants of the molecule

rotational_symmetry_number(self: psi4.core.Molecule) → int

Returns number of unique orientations of the rigid molecule that only interchange identical atoms

rotor_type(self: psi4.core.Molecule) → str

Returns rotor type, e.g. ‘RT_ATOM’ or ‘RT_SYMMETRIC_TOP’

run_dftd3(func=None, dashlvl=None, dashparam=None, dertype=None, verbose=False)

Compute dispersion correction using Grimme’s DFTD3 executable.

Function to call Grimme’s dftd3 program to compute the -D correction of level dashlvl using parameters for the functional func. dashparam can supply a full set of dispersion parameters in the absence of func or individual overrides in the presence of func.

The DFTD3 executable must be independently compiled and found in PATH or PSIPATH.

Parameters:
  • mol (qcdb.Molecule or psi4.core.Molecule or str) – Molecule on which to run dispersion calculation. Both qcdb and psi4.core Molecule classes have been extended by this method, so either allowed. Alternately, a string that can be instantiated into a qcdb.Molecule.
  • func (str or None) – Density functional (Psi4, not Turbomole, names) for which to load parameters from dashcoeff[dashlvl][func]. This is not passed to DFTD3 and thus may be a dummy or None. Any or all parameters initialized can be overwritten via dashparam.
  • dashlvl ({'d2p4', 'd2gr', 'd3zero', 'd3bj', 'd3mzero', d3mbj', 'd', 'd2', 'd3', 'd3m'}) – Flavor of a posteriori dispersion correction for which to load parameters and call procedure in DFTD3. Must be a keys in dashcoeff dict (or a key in dashalias that resolves to one).
  • dashparam (dict, optional) – Dictionary of the same keys as dashcoeff[dashlvl] used to override any or all values initialized by dashcoeff[dashlvl][func].
  • dertype ({None, 0, 'none', 'energy', 1, 'first', 'gradient'}, optional) – Maximum derivative level at which to run DFTD3. For large mol, energy-only calculations can be significantly more efficient. Also controls return values, see below.
  • verbose (bool, optional) – When True, additionally include DFTD3 output in output.
Returns:

  • energy (float, optional) – When dertype is 0, energy [Eh].
  • gradient (list of lists of floats or psi4.core.Matrix, optional) – When dertype is 1, (nat, 3) gradient [Eh/a0].
  • (energy, gradient) (float and list of lists of floats or psi4.core.Matrix, optional) – When dertype is unspecified, both energy [Eh] and (nat, 3) gradient [Eh/a0].

Notes

research site: https://www.chemie.uni-bonn.de/pctc/mulliken-center/software/dft-d3 Psi4 mode: When psi4 the python module is importable at import qcdb

time, Psi4 mode is activated, with the following alterations: * output goes to output file * gradient returned as psi4.core.Matrix, not list o’lists * scratch is written to randomly named subdirectory of psi scratch * psivar “DISPERSION CORRECTION ENERGY” is set * verbose triggered when PRINT keywork of SCF module >=3
run_gcp(func=None, dertype=None, verbose=False)

Function to call Grimme’s dftd3 program (http://toc.uni-muenster.de/DFTD3/) to compute the -D correction of level dashlvl using parameters for the functional func. The dictionary dashparam can be used to supply a full set of dispersion parameters in the absense of func or to supply individual overrides in the presence of func. Returns energy if dertype is 0, gradient if dertype is 1, else tuple of energy and gradient if dertype unspecified. The dftd3 executable must be independently compiled and found in PATH or PSIPATH. self may be either a qcdb.Molecule (sensibly) or a psi4.Molecule (works b/c psi4.Molecule has been extended by this method py-side and only public interface fns used) or a string that can be instantiated into a qcdb.Molecule.

save_string_xyz(self: psi4.core.Molecule) → str

Saves the string of an XYZ file to arg2

save_string_xyz_file(self: psi4.core.Molecule) → str

Saves an XYZ file to arg2

save_xyz_file(self: psi4.core.Molecule, arg0: str, arg1: bool) → None

Saves an XYZ file to arg2

schoenflies_symbol(self: psi4.core.Molecule) → str

Returns the Schoenflies symbol

scramble(do_shift=True, do_rotate=True, do_resort=True, deflection=1.0, do_mirror=False, do_plot=False, run_to_completion=False, run_resorting=False, verbose=1)

Tester for B787 by shifting, rotating, and atom shuffling ref_mol and checking that the aligner returns the opposite transformation.

Parameters:
  • ref_mol (qcdb.Molecule or psi4.core.Molecule) – Molecule to perturb.
  • do_shift (bool or array-like, optional) – Whether to generate a random atom shift on interval [-3, 3) in each dimension (True) or leave at current origin. To shift by a specified vector, supply a 3-element list.
  • do_rotate (bool or array-like, optional) – Whether to generate a random 3D rotation according to algorithm of Arvo. To rotate by a specified matrix, supply a 9-element list of lists.
  • do_resort (bool or array-like, optional) – Whether to shuffle atoms (True) or leave 1st atom 1st, etc. (False). To specify shuffle, supply a nat-element list of indices.
  • deflection (float, optional) – If do_rotate, how random a rotation: 0.0 is no change, 0.1 is small perturbation, 1.0 is completely random.
  • do_mirror (bool, optional) – Whether to construct the mirror image structure by inverting y-axis.
  • do_plot (bool, optional) – Pops up a mpl plot showing before, after, and ref geometries.
  • run_to_completion (bool, optional) – By construction, scrambled systems are fully alignable (final RMSD=0). Even so, True turns off the mechanism to stop when RMSD reaches zero and instead proceed to worst possible time.
  • run_resorting (bool, optional) – Even if atoms not shuffled, test the resorting machinery.
  • verbose (int, optional) – Print level.
Returns:

Return type:

None

set_active_fragment(self: psi4.core.Molecule, arg0: int) → None

Sets the specified fragment arg2 to be Real

set_active_fragments(self: psi4.core.Molecule, arg0: List[int]) → None

Sets the specified list arg2 of fragments to be Real

set_basis_all_atoms(self: psi4.core.Molecule, arg0: str, arg1: str) → None

Sets basis set arg2 to all atoms

set_basis_by_label(self: psi4.core.Molecule, arg0: str, arg1: str, arg2: str) → None

Sets basis set arg3 to all atoms with label (e.g., H4) arg2

set_basis_by_symbol(self: psi4.core.Molecule, arg0: str, arg1: str, arg2: str) → None

Sets basis set arg3 to all atoms with symbol (e.g., H) arg2

set_geometry(self: psi4.core.Molecule, arg0: psi4.core.Matrix) → None

Sets the geometry, given a (Natom X 3) matrix arg2 of coordinates (in Bohr)

set_ghost_fragment(self: psi4.core.Molecule, arg0: int) → None

Sets the specified fragment arg2 to be Ghost

set_ghost_fragments(self: psi4.core.Molecule, arg0: List[int]) → None

Sets the specified list arg2 of fragments to be Ghost

set_input_units_to_au(self: psi4.core.Molecule, arg0: float) → None

Sets unit conversion to [a0] for geometry

set_mass(self: psi4.core.Molecule, atom: int, mass: float) → None

Sets mass of atom (0-indexed) to mass

set_molecular_charge(self: psi4.core.Molecule, arg0: int) → None

Sets the molecular charge

set_multiplicity(self: psi4.core.Molecule, arg0: int) → None

Sets the multiplicity (defined as 2Ms + 1)

set_name(self: psi4.core.Molecule, arg0: str) → None

Sets molecule name

set_nuclear_charge(self: psi4.core.Molecule, arg0: int, arg1: float) → None

Set the nuclear charge of the given atom to the value provided.

set_point_group(self: psi4.core.Molecule, arg0: psi4.core.PointGroup) → None

Sets the molecular point group to the point group object arg2

set_units(self: psi4.core.Molecule, arg0: psi4.core.GeometryUnits) → None

Sets units (Angstrom or Bohr) used to define the geometry

set_variable(self: psi4.core.Molecule, arg0: str, arg1: float) → None

Assigns the value arg3 to the variable arg2 in the list of geometry variables, then calls update_geometry()

symbol(self: psi4.core.Molecule, arg0: int) → str

Gets the cleaned up label of atom arg2 (C2 => C, H4 = H)

symmetrize(self: psi4.core.Molecule, arg0: float) → None

Finds the highest point Abelian point group within the specified tolerance, and forces the geometry to have that symmetry.

symmetry_from_input(self: psi4.core.Molecule) → str

Returns the symmetry specified in the input

to_arrays()

Exports coordinate info into NumPy arrays.

Returns:
  • geom, mass, elem, elez, uniq (ndarray, ndarray, ndarray, ndarray, ndarray) – (nat, 3) geometry [a0]. (nat,) mass [u]. (nat,) element symbol. (nat,) atomic number. (nat,) hash of element symbol and mass. Note that coordinate, orientation, and element information is preserved but fragmentation, chgmult, and dummy/ghost is lost.
  • Usage
  • —–
  • geom, mass, elem, elez, uniq = molinstance.to_arrays()
to_dict(force_c1=False, force_units=False, np_out=True)

Serializes instance into Molecule dictionary.

to_schema(dtype, units='Angstrom', return_type='json')

Serializes instance into JSON or YAML according to schema dtype.

to_string(dtype, units='Angstrom', atom_format=None, ghost_format=None, width=17, prec=12)

Format a string representation of QM molecule.

translate(self: psi4.core.Molecule, arg0: psi4.core.Vector3) → None

Translates molecule by arg2

units(self: psi4.core.Molecule) → str

Returns units used to define the geometry, i.e. ‘Angstrom’ or ‘Bohr’

update_geometry(self: psi4.core.Molecule) → None

Reevaluates the geometry with current variable values, orientation directives, etc. Must be called after initial Molecule definition by string.

x(self: psi4.core.Molecule, arg0: int) → float

x position of atom

y(self: psi4.core.Molecule, arg0: int) → float

y position of atom

z(self: psi4.core.Molecule, arg0: int) → float

z position of atom