Dear openbabel support team,
during the attempt to search for the minimum energy rotational conformer of some nitrogenated and oxygenated species, like they are important for detailed modeling of combustion situations, I'm occuring problems to find correct results with the current set of functions I'm aware of within the python frontend of openbabel. The tool I'm working on shall help to identify rotational conformers as a preceding step towards a higher level quantum mechanical calculation. The focus is set on intermediate to longer chain hydrocarbons with attached functional groups. Those might be nitrogenated species with integrated or added nitrogen, like nitrites, nitrates or hydrazines, and also oxygenated with peroxy or hydroperoxy-groups attached. When it comes to e.g. n-alkylnitrites it starts to get difficult to identify the minimum rotational conformer for further subsequent calculations because of the large degree of rotational freedom. Such a tool would be very helpful for those situations. Within the development process of the tool I initially followed the suggestions of the openbabel manual on single conformer searching like it is described here: https://open-babel.readthedocs.io/en/latest/3DStructureGen/SingleConformer.html Unfortunately, on some test cases, e.g. the ROOHs of acetaldehyde (O=CCOO, O=C(OO)C) the results of two different force fields (gaff, mmff94) were inconsistent and different in the regard of internal hydrogen bonds between the carbonyl group and the hydrogen atom of the hydroperoxy group. Similarily, the two different force fields lead to different results for n-pentylnitrite (Within the gen3d approach). Also for hydrazine (NN) the method is not capable to identify both and more importantly the minimum energy conformer at the internal bond rotation. Just the anti conformer is identified and not the gauche conformer, where the hydrogen bond is formed between the lone-pairs of the nitrogen atoms with one of the opposing hydrogens. Generally, all the available force fields in openbabel and in particular gaff and mmff94 are capable to describe the hydrogen bonded structures satisfyingly of species containing H, C, N and O. Therefore, the reason for the inconsistent results seems to be more statistically within the nature of the weighted rotor search method, and actually also the other ones which I tried as well. The common denominator of all these methods on first glance seems to be the torsional database (torlib.txt), where seemingly no definitions for ...N-N..., ...O-N..., ...O-O... bonds are available. So I tried to alter it (without any knowledge) by adding the following entries: # hydrazines [C,H]N([C,H])~N([C,H])[C,H] 1 2 4 5 180.0 60.0 -60.0 120.0 -120.0 0.0 180.0 30.0 -30.0 # Peroxy, hydroperoxy, nitrites *O[O,N]* 1 2 3 4 180.0 60.0 -60.0 120.0 -120.0 0.0 180.0 30.0 -30.0 But this did not change anything on the results. Maybe you could give me some more tipps or a better description of the smarts and groups in that regard. Because on my lacking knowledge on all of the details of openbabel, I thought it is a good idea to start with a new approach, a brute force analysis of all rotational conformers of the internal bonds, loosely following that example: https://gist.github.com/ghutchis/0d89de9f81157cf1aa9726018d41e8e2 And then I recognized, that not all internal rotational bonds are identified, pretty similar like it was previously also stated from Scott McKechnie https://sourceforge.net/p/openbabel/mailman/message/32125004/ In my case, with the following minimal code snippet: import openbabel #smi = 'O=CCOO' smi = 'NN' # . openbabel objects mol = openbabel.OBMol() conv = openbabel.OBConversion() build = openbabel.OBBuilder() ff = openbabel.OBForceField.FindForceField('gaff') # . initialize molecule conv.SetInAndOutFormats('smi', 'mdl') conv.ReadString(mol, smi) mol.AddHydrogens() build.Build(mol) ff.Setup(mol) ff.SteepestDescent(250) ff.GetCoordinates(mol) # set up the pool of rotatable bonds rl = openbabel.OBRotorList() rl.Setup(mol) print("Number of rotatable bonds: {}".format(rl.Size())) the result for NN gives 0 out of 1 internal rotors and for O=CCOO 2 out of 3 internal rotors. Seemingly, openbabel does not recognize the ...O-O... and the ...N-N... coordinate as an internal bond. Additionally I'm presuming, that the same is happening with ...N-O... bonds, like in the ...ONO or ...ONO2 groups of nitrites and nitrates. For my research it is quite important to include those kind of rotations in the conformational analysis. My suspicion would be, that the inclusion of such rotors also solves the problems on the different results with different forcefields for these cases. Also for rotational scans with openbabel, identifying single rotational bonds in small molecules, like NN, or CC could be helpful. Could you please give me suggestions, hints and ideas to persuade openbabel to consider this kind of structures/rotations? Best regards, Heiko Minwegen
_______________________________________________ OpenBabel-discuss mailing list OpenBabel-discuss@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/openbabel-discuss