Computational Chemistry 1.9 - Torsion Angles

Short lecture on torsion angles in molecules.
If in 4 atoms atoms 1 and 2, 2 and 3, and 3 and 4 are bonded to one another, then these 4 atoms form a torsion angle. The torsion angle is the angle between the normal vectors of plane 123 and plane 234.
Notes Slide: i.imgur.com/Qpz...
-- About TMP Chem --
All TMP Chem content is free for everyone, everywhere, and created independently by Trent Parker.
Email: [e-mail]
-- Video Links --
Course Playlist: • Computational Chemistry
Other Courses: • PChem Course Intros
Channel Info: • About TMP Chem
-- Social Links --
Facebook: / tmpchem
Twitter: / tmpchem
LinkedIn: / tmpchem
Imgur: tmpchem.imgur.com
GitHub: www.github.com...
-- Equipment --
Microphone: Blue Yeti USB Microphone
Drawing Tablet: Wacom Intuos Pen and Touch Small
Drawing Program: Autodesk Sketchbook Express
Screen Capture: Corel Visual Studio Pro X8

[@jonting776]

There was a mistake in the formula to find the normal to the plane abc, cross instead of dot. But it is correct when you use it in your example.

[@TMPChem]

Correct. Thanks for pointing that out. I was very disappointed when I clicked on this image and saw that it actually was incorrect. Oof.

[@jonting776]

No prob. Thanks for the hard works tho! Really beneficial to newbies to computational chemistry like me.:)

[@rishabghosh7685]

Thank you very much for these videos. You are a life saver.

[@andresfelipecanobotero8814]

Hi! First of all, thank you for your work! It is helping me a lot getting started with CC! I wanted to ask if there is a reason why the sin(theta_ijk)*sin(theta_jkl) is not dividing the dot product of r's in your code (what you defined as 'dp')

[@Nitroxien]

Hey Andres, responding since unfortunately TMP Chem has been gone for a while now (hoping he will prove me wrong and comment and say he is coming back soon :) ) however I think your getting stuck on the mistake he made where here defined the normal vector to the plane abc using dot product instead of cross product. The equation should be rba x rbc/sin(theta) = n

[@faycalaissabrahim8036]

i always get a frame of 1 after minimizing in avogadro !!! but in this video its 636 !!! so i get no dynamics !!! help me plz !!!

[@robertpapoular9612]

Happy New Year 2021 ! ... and many thanks for your intriguing video. I would most humbly suggest that ... N_ijk . r_ji ... is most likely to be ... exactly ... zero, since N_ijk is the cross product of r_ji and r_jk. ... or I am getting it all wrong ? If so, where do we get the sign of PHI_ijkl ?

[@TMPChem]

Hi Robert. You are correct. That is a typo. It should be r_kl. Thanks for noting the correction.

[@robertpapoular9612]

@@TMPChem Thank you very much for your feedback. Everything is clear to me now. Best,

[@Greenbeanandcheese]

Not to be nitpicky, but for the cross product equation, shouldn't it be -(azbx-axbz)y instead of +(azbx-axbz)?

[@TMPChem]

Hi Eduardo. I believe the expression is correct as written. There are two different widely-used methods for generating the terms in a cross product. The first is to do determinant expansion, in which case as you translate to the right, you would alternate between + and - as you noted, where you only include the columns which aren't part of the current unit vector. The second is to do a cyclic permutation of all indices (x -> y, y -> z, z -> x), in which case you the sign doesn't alternate, you just copy the result for the previous dimension and permute the indices, repeating until the indices have wrapped around.

[@luosong4421]

I'm sorry but why the torsion angle of HOOH was -180°, it seems like it was quite different with what we saw in VMD

[@TMPChem]

Hi Luo. You are correct the coordinates that are being analyzed in the hooh.xyz file from the geom/xyz directory in the GitHub repository are different than the ones used to generate that trajectory image in VMD. The hydrogen peroxide coordinates analyzed by the torsion.py program had hydrogen atoms that are trans relative to one another, giving a torsion angle of -180.0 degrees.

[@luosong4421]

Thank you for your response. Now, I get it.Thanks a lot

[@1234skyshark]

What can reduce the atom quartet comparisons to O(n log n)?

[@TMPChem]

Appropriate use of algorithms and data structures, which can get quite involved and abstract, but I'll try to summarize (not-so) briefly. The first key is to note that only cases in which atom pairs 12, 23, and 34 are bonded are considered to be torsions. So if 1 and 2 aren't bonded, we need not compare any quartets which necessitate 1 and 2 being bonded (and analogously for 23 and 34). This will already on its own reduce the number of comparisons to O(n^2). Note that in the program I wrote, the first step in determining the torsions is to build a bond graph, which is an array / list of all atom pairs which are bonded. Naively it costs O(n^2) to build this array. Each individual atom has a maximum number of bonds, which is a constant [ O(1) ], and there are n of them. Thus, once I have the bond graph it costs O(n) to search it for torsions as is done in the get_torsions function. Yes, it has a quadruply nested for loop, but the three inner loops have a constant number of elements, so only the outer scales as O(n). Since it costs O(n^2) to build the bond graph, this gives us an overall scaling of O(n^2) to find all of the valid torsions among all O(n^4) atom quartets. So the question of reducing torsion discovery to O(n log n) is really equivalent to asking how we can build the bond graph using an algorithm which scales as O(n log n) instead of doing O(n^2) atom pair comparisons. This is where things get tricky. First, let's note that if (x2 - x1), (y2 - y1), *or* (z2 - z1) are greater than the threshold distance [ k (r_2 + r_1), see the covalent radius video ], then the two atoms are not going to be bonded. Next, let's find the atom in the system with the *largest* covalent radius, and note that this is still true if instead of [ k (r_2 + r_1) ], we use [k (r_max + r_max) = 2 k r_max ]. Typically r_max is < 2 Angstroms, and 2 k r_max is less than 3 Angstroms, so let's use 3 Angstroms as our "key distance". Next, let's partition all of space into blocks of 3 x 3 x 3 Angstroms, and assign every atom to one of these boxes. We'll then assign a key to every atom based on which block it's in. We then build an array which contains a list of every atom which is contained inside every block. Since atoms have a finite size, there is a maximum number of atoms which can be inside each block, typically around 5. Next, we note that any block which is not a direct neighbor to the block our atom is in cannot be bonded to our atom. So we only need to compare atoms that are in neighboring blocks. There are 26 total neighbors, 6 which share a face, 12 which share an edge, and 8 which share a corner. Each atom belongs to one block with