Local force constants

Local force constants can be calculated with the local modes method [1] or the compliance matrix method [2]. ᴍᴏʀғᴇᴜs can use the output of the Gaussian, xtb, or UniMoVib programs.

Preparing input

The LocalForce class needs input from quantum-chemical frequency calculations. Therefore, the instructions are a bit more involved. This input can be read using the load_file method, and ᴍᴏʀғᴇᴜs can also calculate some quantities using its built-in normal mode analysis and internal coordinate codes.

Local modes

The local modes approach needs the normal modes in internal coordinates, here called internal modes, as well as the normal mode force constants. The internal modes can either be read from file or calculated with ᴍᴏʀғᴇᴜs. Specifically, the internal modes can be computed from the normal modes and the Wilson B matrix. For local frequencies, additional input is needed.

Input needed for the local modes approach.

Compliance matrix

The compliance matrix method needs the Hessian matrix and the Wilson B matrix. For local frequencies, the elements are also needed (for their atomic masses).

Input needed for the compliance matrix method.

Recommended input

The table below summarizes the recommended input to ᴍᴏʀғᴇᴜs from the supported programs. More details are given below.

Recommended input.

Method	Gaussian	UniMoVib	xtb
Local modes	log	log	hessian
Compliance	fchk	umv/local	hessian

Gaussian

Gaussian can provide the complete input for both the local modes and the compliance matrix methods. The log file is recommended as input for the the local modes method, and the fchk file for the compliance matrix method. The Gaussian input file should be prepared with the following keywords included (appropriate method and basis set need to be added as well).

log	freq=intmodes iop(7/75=-1) iop(1/33=3)
fchk	freq=intmodes

iop(7/75=1) triggers printing of the internal modes with full accuracy, while iop(1/33=3) triggers printing of the Wilson B matrix. The table below lists all the quantities read by each approach. fchk files should be generated with formchk -3 to obtain maximum printout of internal coordinates.

Quantities from Gaussian files.

Quantity	log	fck
Elements	x	x
Coordinates	x	x
Internal modes	x
Internal coordinates	x	x
B Matrix	x
Normal modes		x
Hessian		x
Force constants	x	x

Note

Additional accuracy for local modes calculations can be achieved by loading also the fchk file which contains the normal mode force constants to higher accuracy.

UniMoVib

The UniMoVib program can do vibrational analysis for a range of different quantum-chemical programs. The user needs to define internal coordinates in addition to the information that can be read from the files. For the local modes method, the log file is recommended. For the compliance matrix method, either the umv or the local files are recommended.

Quantities from UniMoVib files.

Quantity	local	log	umv
Elements	x	x	x
Coordinates	x	x	x
Internal modes
Internal coordinates
B Matrix
Normal modes	x	x
Hessian	x		x
Force constants		x

The log file is the standard output of the program. The umw file and the local file can be generated by specifying them in the input file:

Example UniMoVib input file

a test job

$contrl
  qcprog="gaussian"
  iflocal=.t.
  ifsave=.t.
$end

$qcdata
  fchk="freq.fchk"
$end

xtb

The xtb program can provide the Hessian. The file hessian is generated by the xtb program as a result of a frequency calculation. It needs to be used together with a normal mode analysis and separate internal coordinate input.

Quantities from xtb files.

Quantity	hessian
Elements
Coordinates
Internal modes
Internal coordinates
B Matrix
Normal modes
Hessian	x
Force constants

Geometry optimizations

Local force constants and frequencies are based on the harmonic approximation that is valid only at stationary points (minima and transition states) on the potential energy surface. Therefore, the geometry optimizations need to be of good quality. For difficult examples on flat surfaces, this can mean increasing the convergence criteria of the quantum-chemical program. Vibrations with small imaginary frequencies should be eliminated as much as possible. The local modes method with standards settings is somewhat robust to the presence of these vibrations , while they can result in artifacts with the compliance matrix method. Therefore, the local modes method is recommended in these cases.

For transition states, the imaginary mode corresponding to the reaction is projected out with the local modes approach. This also means that local force constant involving the atoms corresponding to this imaginary mode are meaningless and should not be used. Only the local modes method can be used for transition states.

Warning

Transition states can only be treated with the local modes method. Force constants and frequencies of the atoms corresponding to the imaginary mode should not be used.

Module

The LocalForce class is provided to calculate and store the local force constants and frequencies.

Example with Gaussian and local modes method.

>>> from morfeus import LocalForce
>>> lf = LocalForce()
>>> lf.load_file("freq-lm.log", "gaussian", "log")
>>> lf.compute_local()
>>> lf.compute_frequencies()
>>> fc = lf.get_local_force_constant([1, 2])
>>> print(fc)
5.364289643211871
>>> freq = lf.get_local_frequency([1, 2])
>>> print(freq)
3129.3126301763527

Example with Gaussian and compliance matrix method.

>>> from morfeus import LocalForce
>>> lf = LocalForce()
>>> lf.load_file("freq.fchk", "gaussian", "fchk")
>>> lf.compute_compliance()
>>> lf.compute_frequencies()
>>> fc = lf.get_local_force_constant([1, 2])
>>> print(fc)
5.364398642985929
>>> freq = lf.get_local_frequency([1, 2])
>>> print(freq)
3129.352986019491

Example with xtb and local modes method.

>>> from morfeus import LocalForce, read_xyz
>>> elements, coordinates = read_xyz("xtbopt.xyz")
>>> lf = LocalForce(elements, coordinates)
>>> lf.load_file("hessian", "xtb", "hessian")
>>> lf.normal_mode_analysis()
>>> lf.detect_bonds()
>>> print(lf.internal_coordinates)
[Bond(1, 4), Bond(1, 3), Bond(1, 2), Bond(1, 5)]
>>> lf.compute_local()
>>> lf.compute_frequencies()
>>> fc = lf.get_local_force_constant([1, 2])
>>> print(fc)
5.190222259808879
>>> freq = lf.get_local_frequency([1, 2])
>>> print(freq)
3078.130379468432

Example with UniMoVib and the local modes method.

>>> from morfeus import LocalForce
>>> lf = LocalForce()
>>> lf.load_file("job.out", "unimovib", "log")
>>> lf.detect_bonds()
>>> lf.compute_local()
>>> lf.compute_frequencies()
>>> fc = lf.get_local_force_constant([1, 2])
>>> print(fc)
5.364347084281302
>>> freq = lf.get_local_frequency([1, 2])
>>> print(freq)
3129.337947449028

Example with adding internal coordinates manually

>>> from morfeus import LocalForce
>>> lf = LocalForce()
>>> lf.load_file("job.out", "unimovib", "log")
>>> lf.add_internal_coordinate([1, 2])
>>> lf.add_internal_coordinate([1, 2, 3])
>>> print(lf.internal_coordinates)
[Bond(1, 2), Angle(1, 2, 3)]
>>> lf.compute_local()
>>> lf.compute_frequencies()
>>> fc = lf.get_local_force_constant([1, 2])
>>> print(fc)
5.364347084281298
>>> freq = lf.get_local_frequency([1, 2])
>>> print(freq)
3129.337947449028
>>> lf.print_report(angles=True, angle_units=True)
Coordinate                            Force constant (mDyne/Å, mDyne Å rad^(-2))             Frequency (cm^-1)
Bond(1, 2)                                                                 5.364                          3129
Angle(1, 2, 3)                                                             2.416                          1687

For the local modes method, projection of imaginary frequencies can be controlled with the project_imag=<bool> keyword to the compute_local method. The cutoff for low-frequency modes can be controlled with cutoff=<float>. Internal coordinates can be added with the add_internal_coordinate method.

For more detailed information, use help(LocalForce) or see the API documentation: LocalForce.

Command line script

The command line script provides access to the basic functionality through the terminal.

Example with Gaussian log file.

$ morfeus local_force - - load_file freq-lm.log gaussian log - compute_local - compute_frequencies - print_report
Coordinate                                              Force constant (mDyne/Å)              Frequency (cm⁻¹)
Bond(1, 2)                                                                 5.364                          3129
Bond(1, 3)                                                                 5.364                          3129
Bond(1, 4)                                                                 5.365                          3129
Bond(1, 5)                                                                 5.364                          3129

Note

Breaking up the command:

morfeus local_force - -

Call LocalForce without any geometry file and without arguments. Return object.

load_file freq_lm.log gaussian log -

Call LocalForce.load_file with three arguments. Return object for chaining.

compute_local -

Call LocalForce.compute_local. Return object for chaining.

compute_frequencies -

Call LocalForce.compute_frequencies. Return object for chaining.

print_report

Call LocalFroce.print_report. Return result

Another example using xtb where a geometry file is loaded:

Example with xtb and hessian file.

$ morfeus local_force xtb.xyz - - load_file hessian xtb hessian - normal_mode_analysis -  detect_bonds - compute_local - compute_frequencies - print_report
Coordinate                                              Force constant (mDyne/Å)              Frequency (cm⁻¹)
Bond(1, 2)                                                                 5.190                          3078
Bond(1, 3)                                                                 5.190                          3078
Bond(1, 4)                                                                 5.190                          3078
Bond(1, 5)                                                                 5.190                          3078

Background

Local force constants describe the bond strength based on vibrational frequencies. There are two approaches in the literature: the local modes method of Cremer [1] and the compliance matrix method of Grunenberg [2]. They have been shown to be equivalent within numerical accuracy [3]. ᴍᴏʀғᴇᴜs can use either method, and they give almost identical results except when there are modes with imaginary or very small frequencies. In these cases, the local modes approach can handle the issue with two methods: (1) projecting out imaginary modes, and (2) raising the force constants of low-frequency modes to a cutoff value. ᴍᴏʀғᴇᴜs does this projection by default and uses a cutoff of 0.001 mDyne/Å for low-frequency modes. We therefore recommend local modes with default settings as the most robust method in problematic cases. Expert users can turn off the projection and alter the cutoff value. Note that interactions involving imaginary modes (such as breaking/forming bonds in transition states) cannot be assessed by the local force constants. ᴍᴏʀғᴇᴜs has been benchmarked against the local force constants and frequencies given by Cremer [3] for small organic molecules.

References

[1] (1,2)

Zoran Konkoli and Dieter Cremer. A new way of analyzing vibrational spectra. I. Derivation of adiabatic internal modes. International Journal of Quantum Chemistry, 67(1):1–9, 1998. doi:10.1002/(SICI)1097-461X(1998)67:1<1::AID-QUA1>3.0.CO;2-Z.

[2] (1,2)

Kai Brandhorst and Jörg Grunenberg. How strong is it? The interpretation of force and compliance constants as bond strength descriptors. Chemical Society Reviews, 37(8):1558–1567, 2008. doi:10.1039/B717781J.

[3] (1,2)

Wenli Zou, Robert Kalescky, Elfi Kraka, and Dieter Cremer. Relating normal vibrational modes to local vibrational modes with the help of an adiabatic connection scheme. The Journal of Chemical Physics, 137(8):84114, 2012. doi:10.1063/1.4747339.