Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Feedback on SRS #34

Open
PaulWAyers opened this issue Nov 29, 2018 · 0 comments
Open

Feedback on SRS #34

PaulWAyers opened this issue Nov 29, 2018 · 0 comments

Comments

@PaulWAyers
Copy link

  1. The Global minimum is almost always meaningless because it is not a conformational isomer but a true isomer (where bonds are broken). So we want to restrict our conformational search to rotations around bonds. (One could go further, including things like H-atom transfers across hydrogen bonds and other similar low-barrier processes, but that would make things significantly more complicated.)

  2. You can be a bit more clear about what you are doing. You are searching for the best choice of bond torsions/dihedral angles. (Of course there will be small changes in bond angles/lengths associated with these changes, but they are induced by the torsions.) Realistically, you may wish to construct the whole ensemble (all the thermally accessible structures) not just the lowest-energy conformer. That may be beyond the mandate of the program, but need not be....usually you will find many (if not all) of the other low-energy conformers in the process of seeking the lowest-energy conformer.

  3. Is it essential that the user specify (guess) the number of conformers? Can this be guessed by some heuristic? I would think that the number of rotatable bonds could be identified, and then multiplied by ~3, to get an upper bound on the number of conformers. (The number of low-energy conformers will, in general, be much smaller than this....something of the order of N_atoms or N_atoms^2 I would guess, except (again always an exception) in painfully "glassy" molecules (e.g., long floppy polymers) where the problem is more-or-less intractable (and uninteresting as many conformers are chemically equivalent) also.

  4. 5.1 The RMSD is tricky. Sometimes there are, due to symmetry, several equivalent conformers. I am not sure it is practical (it might be very hard) but technically while only one of these conformers needs to be stored, the "degeneracy" (number of equivalent instances of this conformer) should be also stored, so that one can compute thermodynamic properties. Again, this is if the mandate is to allow thermal averages over conformers (in which case the degeneracy of each conformer is needed) to compute thermodynamic quantities. If the goal is merely to identify the lowest-energy conformer(s) for subsequent studies (qualitative reactivity and/or quantitative electronic/mechanical properties studies) then this is not essential. It isn't clear to me exactly what the mandate is....just find some low-energy conformers that are inequivalent? Is the temperature just a screening criterion or is it intended to be used to guide thermal averaging later? The Procrustes package (Fanwang) has RMSD and a lot of other things like that already there, though the RMSD between molecules is pretty straightforward as long as you don't worry about matching atom types (and only atomic positions).

  5. Conformational isomer. Rotations around triple and quintuple bonds are also free; basically all rotations around bonds where there are an even number of occupied pi and/or delta orbitals. That mostly happens in transition-metal compounds, though, as conformations around triple (or higher) bonds require atoms that form at least 5 bonds. Still, what you say isn't quite precisely correct... It is probably better to merely say that single bonds are rotatable and then the user can specify (by hand) other bonds they wish to consider rotatable, because there are cases where there are (weak) double bonds that are also rotatable.

  6. I'm not sure I trust A8. Usually you need to converge the bond lengths/angles at least loosely. You may not need to do this in the very early stages but, for example, the rotational barrier in simple hydrocarbons (especially if there are bulky-ish substituents) is quite far off if you freeze bond lengths and angles. Something similar would be true of chair/boat/twisted-boat/etc. geometries for ring compounds. Probably it is safe to do initial screening without optimization and after that a coarse optimization with early stopping once you have converged the geometry to an energy that is small enough (say 1/2 of the conformer energy gap you are trying to assess) is adequate.

  7. Page 9. The nuclei are much more massive than the electrons. They aren't much larger...if fact they are about 100,000 times smaller (in terms of the spatial extent of their wavefunctions).

  8. The AOs are pretty muddled. You normally write something like cexp(-ar^{2})*r^(some power)*spherical harmonic. (Obviously not r^{2} in the exponential for Slater orbitals.) Slater orbitals do not diverge at the nucleus (no -> infinity) but they have a derivative discontinuity/cusp at the nucleus, which is "correct" if we assume that the nucleus is a point charge. You are right that it is easier to do the integrals for Gaussian orbitals and that while the near-nucleus and far-from-nucleus behavior of Gaussians is generally poorer than for Slaters, taking a linear combination of Gaussians to approximate a Slater type function is a pragmatic compromise. In practice, doing the integrals should not (but this is not always true) be the bottleneck in your code anyway.....

  9. Page 11. The RMSD cannot be computed until the optimal alignment (align centers of mass of the two conformers and then rotate the conformers to minimize the RMSD; that is the thing that Fanwang's code does though again you can safely implement it in situ) is constructed. So the formula you give is a RMSD but not what we normally mean by the RMSD between two molecules/conformers.

  10. During the process of the exploration S_{RMSD} should be reduced. Early in the optimization it is favorable to force yourself to explore very diverse structures. Later on you want to take the best...even if they are not very far from each other. Alternatively make the RMSD criterion nonlinear, something like (1-exp(-RMSD^2)) so that as soon as two structures are "different enough" they contributed a "big enough" amount to the diversity criterion.

  11. Appendix. Misspelled "Planck" and also should maybe (not clear to me from the text what you intend here) include the nuclear-nuclear repulsion energy in the Hamiltonian.

Overall: Looks very good. As a specification I like it, though I'm not sure algorithmically....achieving the specification may be tricky.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant