0
$\begingroup$

Is there any software that I can provide bond lengths, bond angles and torsion dihedrals of a molecule and generate an output file in the format of cartesian coordinates or zmatrix?

Example:

I would have input similar to the one below (which i provide all bond lengths, bond angles and torsion dihedrals of a molecule)

BOND LENGTHS
     2 1 Si O2 1.63
     3 1 Si O2 1.63
     ...
BOND ANGLES
     2 1 3 Si O2 Si 148,482
     1 2 4 O2 Si C3 107,426
     ...
TORSION ANGLES
     3 1 2 4 -161,552
     3 1 2 8 -41,721
     ...

And I would like to have an output like this:

Example

My idea is to vary some of this information and generate several Gaussian software inputs to calculate energy points at specific bond lengths, for example, without having to manually edit the molecule structure in this software.

I know that Open Babel software can do something like that, but I don't know what type of input to use.

$\endgroup$
  • $\begingroup$ OBabel can write in Fenske-Hall Z-Matrix format, but it cannot read it —going from internal coordinates to cartesian coordinates is very uncommon. It sounds like your starting data would need to be processed anyway. The algorithm (called Natural Extension Reference Frame (NeRF) algorithm) to do the conversion is a few trigonometric functions that would be very straightforward to implement with numpy and rdkit in Python. $\endgroup$ – Matteo Ferla Jul 31 '20 at 10:19
  • $\begingroup$ Thanks for the comment @MatteoFerla! I know OpenBabel, I tried to post a topic in the software group but without success (it seems that it is stopped). My idea is to take a structure of a molecule and vary, for example, the bond length, and thereby calculate some energy points in the Gaussian software. The only way I know how to do that today is to manually change the structure of the molecule in softwares like Avogadro or GaussView. I would like to know if there is any way to generate the inputs of the distorted molecule in the bond without doing this manually by editing its structure. $\endgroup$ – Emerson P L Jul 31 '20 at 11:00
  • $\begingroup$ What software group? The problem with internal coordinates is that they are redundant meaning that there are more internal coordinates than internal degrees of freedom. So there is no 1 to 1 transformation. A z matrix is non redundant however so you can use that to vary your structures. Scans along internal coordinates are frequently performed using z matrices $\endgroup$ – Cody Aldaz Aug 6 '20 at 14:13
1
$\begingroup$

This is not the answer to your question, but it is a comment about an approach to the aims, based on the comments. If coding is an option, RDKit is a good option, specifically adding constraints (restraints) to the ForceField and optimising. Given a molecule:

from rdkit import Chem
from rdkit.Chem import AllChem

acetate = Chem.MolFromSmiles('CC(=O)C')

Normally the 3D conformer is generated with:

AllChem.EmbedMolecule(acetate)

And it stops there. However, then you can manually control the force field —either UFF (more elements) or MMFF (classic choice).

ff = AllChem.UFFGetMoleculeForceField(mol)
ff.Initialize()
print(ff.Minimize()) # this tells you if it failed.

Within this, you can add constraints, such as:

ff.AddDistanceConstraint(0, 1, 2., 2.1, 1_000_000_000)

Note that the weight needs to be set absurdly high for absurd lengths. To check it work you could do:

dm = Chem.Get3DDistanceMatrix(acetate) print(f'{dm[0,1]} Å. {ff.CalcEnergy()} kcal/mol')

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.