Geometry¶
The first step of any quantum chemistry calculation is specification of the system geometry. The InQuanto
GeometryMolecular and GeometryPeriodic classes
exist to standardize geometry formats and provide convenient
functionality for building and manipulating molecular and periodic structures.
Molecular Systems¶
Initializing Structures¶
Geometry objects may be constructed in several ways. Atom-by-atom construction is possible
with the add_atom() method:
from inquanto.geometries import GeometryMolecular
g = GeometryMolecular()
g.add_atom("H", [0, 0, 0])
g.add_atom("H", [0, 0, 0.735])
g.add_atom("H", [0, 0, 1.])
print(g.df)
Element X Y Z
id
0 H 0 0 0.000
1 H 0 0 0.735
2 H 0 0 1.000
where the dataframe attribute is a
pandas.DataFrame of the
geometric structure. If the atom position is left empty in add_atom(),
random coordinates between ±1 Angstroms are generated. See below:
g = GeometryMolecular()
g.add_atom("H")
print(g.df)
Element X Y Z
id
0 H 0.233297 -0.633343 -0.68505
Note
InQuanto Geometry objects support Angstrom (default) and Bohr as distance units, specified by the
distance_units constructor argument. If another unit is used, the user can convert the geometry to Bohrs or Angstroms
using the rescale_position_vectors() class method.
Alternatively, one may instantiate a GeometryMolecular object from z-matrices or
Cartesian coordinates provided they are appropriately formatted. Examples for both are given below.
water_zmatrix = """h
o 1 1.0
h 2 1.0 1 104.5
"""
g = GeometryMolecular(water_zmatrix)
print("Geometry from z-matrix:\n", g.df)
water_xyz = [
["O", [0.0000000, 0.0000000, 0.1271610]],
["H", [0.0000000, 0.7580820, -0.5086420]],
["H", [0.0000000, -0.7580820, -0.5086420]],
]
g = GeometryMolecular(water_xyz)
print("\nGeometry from xyz:\n", g.df)
Geometry from z-matrix:
Element X Y Z
id
0 h 0.0 0.000000 0.00000
1 o 0.0 0.000000 1.00000
2 h 0.0 0.968148 1.25038
Geometry from xyz:
Element X Y Z Atom
id
0 O 0.0 0.000000 0.127161 O1
1 H 0.0 0.758082 -0.508642 H2
2 H 0.0 -0.758082 -0.508642 H3
It is also possible to instantiate/write from files for a limited number of formats without the use of any extensions. These include:
Standard format .XYZ files using the
load_xyz()andsave_xyz()methods..csv files using the
load_csv()andsave_csv()methods..json files using the
load_json()andsave_json()methods.Z-matrix files in the Gaussian format using the
load_zmatrix()andsave_zmatrix()methods.
Calculating and Modifying Properties¶
A variety of geometrical quantities, including bond lengths, angles and dihedrals, can be calculated with
Geometry class methods:
ethane_xyz = [
["C", [0.0000000, 0.0000000, 0.7688350]],
["C", [0.0000000, 0.0000000, -0.7688350]],
["H", [0.0000000, 1.0157000, 1.1533240]],
["H", [-0.8796220, -0.5078500, 1.1533240]],
["H", [0.8796220, -0.5078500, 1.1533240]],
["H", [0.0000000, -1.0157000, -1.1533240]],
["H", [-0.8796220, 0.5078500, -1.1533240]],
["H", [0.8796220, 0.5078500, -1.1533240]],
]
g = GeometryMolecular(ethane_xyz)
print("C=C bond length:", g.bond_length([0, 1]))
print("H-C-C bond angle:", g.bond_angle([2, 0, 1]))
print("2-0-1-5 dihedral angle:", g.dihedral_angle([2, 0, 1, 5]))
g = GeometryMolecular(water_xyz)
print("\n water distance matrix:\n", g.compute_distance_matrix())
C=C bond length: 1.53767
H-C-C bond angle: 110.73395111184878
2-0-1-5 dihedral angle: 180.0
water distance matrix:
[[0. 0.98941082 0.98941082]
[0.98941082 0. 1.516164 ]
[0.98941082 1.516164 0. ]]
These quantities can also be modified:
g = GeometryMolecular(water_xyz)
print("Equilibrium structure 0-1 bond length:", g.bond_length([0, 1]))
print("Equilibrium structure 0-1 bond angle:", g.bond_angle([1, 0, 2]))
# modifying a bond length
g.modify_bond_length(atom_ids=[0, 1], new_bond_length=2)
print("Stretched structure 0-1 bond length:", g.bond_length([0, 1]))
# modifying a bond angle
g.modify_bond_angle(atom_ids=[1, 0, 2], theta=90, units="deg")
print("Closed bond angle 2-0-2:", g.bond_angle([1, 0, 2]))
Equilibrium structure 0-1 bond length: 0.989410821414947
Equilibrium structure 0-1 bond angle: 100.0269116572392
Stretched structure 0-1 bond length: 2.0
Closed bond angle 2-0-2: 90.0
Note
Only the last specified atom in the atom_ids has its position changed when modifying geometric
properties as above. The position of all other atoms remains unchanged.
Calculating and Modifying Properties By Group¶
To transform the geometry by moving more than one atom at a time, one may define a grouping scheme within a structure. A
chemically relevant example is switching between the staggered and eclipsed geometries of ethane. Below we define a
"methyl" grouping scheme which contains two sub-groups:
g = GeometryMolecular(ethane_xyz)
g.set_groups("methyl", {"ch3_1": [0, 2, 3, 4], "ch3_2": [1, 5, 6, 7]})
print(g.df)
Element X Y Z Atom methyl
id
0 C 0.000000 0.00000 0.768835 C1 ch3_1
1 C 0.000000 0.00000 -0.768835 C2 ch3_2
2 H 0.000000 1.01570 1.153324 H3 ch3_1
3 H -0.879622 -0.50785 1.153324 H4 ch3_1
4 H 0.879622 -0.50785 1.153324 H5 ch3_1
5 H 0.000000 -1.01570 -1.153324 H6 ch3_2
6 H -0.879622 0.50785 -1.153324 H7 ch3_2
7 H 0.879622 0.50785 -1.153324 H8 ch3_2
Transformations may then be performed “by group”, so that all atoms within a sub-group are modified simultaneously,
preserving internal structure. For example, below we stretch the C-C bond, then rotate the atoms in "ch3_2"
such that the ethane is in an eclipsed configuration:
g.modify_bond_length_by_group(atom_ids=[0, 1], bond_length=1.8, group="methyl")
g.modify_dihedral_angle_by_group(atom_ids=[2, 0, 1, 7], theta=0.0, group="methyl")
print(g.df)
Element X Y Z Atom methyl
id
0 C 0.000000e+00 0.00000 0.768835 C1 ch3_1
1 C 0.000000e+00 0.00000 -1.031165 C2 ch3_2
2 H 0.000000e+00 1.01570 1.153324 H3 ch3_1
3 H -8.796220e-01 -0.50785 1.153324 H4 ch3_1
4 H 8.796220e-01 -0.50785 1.153324 H5 ch3_1
5 H 8.796220e-01 -0.50785 -1.415654 H6 ch3_2
6 H -8.796220e-01 -0.50785 -1.415654 H7 ch3_2
7 H 3.349405e-16 1.01570 -1.415654 H8 ch3_2
Similarly to single atom modifications, the sub-group that gets moved when transformations are applied is that which contains the last
specified atom in atom_ids.
Generating Many Structures¶
Convenience methods for generating a series of structures with varying geometric properties are also available. For
example, below we generate a range of GeometryMolecular objects for water with
different H-O-H bond angles:
water_geom = GeometryMolecular(water_xyz)
geometries = water_geom.scan_bond_angle([1, 0, 2], [100, 101, 103, 104, 105, 106, 107, 108])
print([g.bond_angle([1,0,2]) for g in geometries])
[np.float64(100.0), np.float64(101.0), np.float64(103.0), np.float64(103.99999999999999), np.float64(105.0), np.float64(106.0), np.float64(107.0), np.float64(108.0)]
Similar functions are implemented for scanning over bond lengths and dihedrals (scan_bond_length()
and scan_dihedral_angle()), and equivalently for varying properties by group
(scan_bond_angle_by_group() for example).
Periodic systems¶
For periodic systems there exists the GeometryPeriodic class, which holds unit cell lattice vectors
in addition to atomic positions:
from inquanto.geometries import GeometryPeriodic
import numpy as np
# Diamond fcc
d = 3.576
a, b, c = [
np.array([0, d/2, d/2]),
np.array([d/2, 0, d/2]),
np.array([d/2, d/2, 0])
]
atoms = [
["C", [0, 0, 0]],
["C", a/4 + b/4 + c/4]
]
cfcc_geom = GeometryPeriodic(geometry=atoms, unit_cell=[a, b, c])
print("Diamond fcc lattice vectors:\n", cfcc_geom.unit_cell)
print("\nCartesian atomic positions:\n", cfcc_geom.df)
Diamond fcc lattice vectors:
[[0. 1.788 1.788]
[1.788 0. 1.788]
[1.788 1.788 0. ]]
Cartesian atomic positions:
Element X Y Z Atom
id
0 C 0.000 0.000 0.000 C1
1 C 0.894 0.894 0.894 C2
GeometryPeriodic supports all of the methods discussed above, excluding instantiation by z-matrices.
In addition, one may easily construct supercells with the build_supercell() method,
which edits the geometry in-place:
cfcc_geom.build_supercell(dimensions=[2, 2, 1])
print("Supercell atomic positions:\n", cfcc_geom.df)
Supercell atomic positions:
Element X Y Z Atom
id
0 C 0.000 0.000 0.000 C1
1 C 0.894 0.894 0.894 C2
2 C 1.788 0.000 1.788 C3
3 C 2.682 0.894 2.682 C4
4 C 0.000 1.788 1.788 C5
5 C 0.894 2.682 2.682 C6
6 C 1.788 1.788 3.576 C7
7 C 2.682 2.682 4.470 C8
To visualize both molecular and periodic structures, see the inquanto-nglview extension.