7.1. Flexibility

In this section, we will describe how to set molecular flexibility for geom.

7.1.1. topgen

There is a small tool called topgen in ABCluster that can be used to generate bonding very conveniently. Assume you have an organic molecule aip.xyz that is in standard XYZ format:

aip.xyz
  1113
  2AIP(+1)
  3N  10.45900000  -12.94900000   -8.94000000
  4H  10.14300000  -12.08200000   -9.34300000
  5H  11.09900000  -12.80100000   -8.17700000
  6H   9.66400000  -13.46100000   -8.55100000
  7C  11.12800000  -13.81000000   -9.92500000
  8H  10.48100000  -13.86000000  -10.80200000
  9C  11.21100000  -15.19900000   -9.28200000
 10O  11.48200000  -15.27400000   -8.08600000
 11C  12.52800000  -13.25600000  -10.32800000
 12H  13.18900000  -13.32500000   -9.46100000
 13C  12.50400000  -11.77400000  -10.78600000
 14H  11.61400000  -11.56800000  -11.38000000
 15H  13.33800000  -11.59800000  -11.46700000
 16C  13.18200000  -14.11500000  -11.43200000
 17H  14.17300000  -13.74000000  -11.68900000
 18H  12.58500000  -14.11600000  -12.34500000
 19H  13.32000000  -15.15200000  -11.12300000
 20N  10.86100000  -16.25800000  -10.02800000
 21H  10.63200000  -16.06600000  -10.99200000
 22C  12.65700000  -10.74000000   -9.65700000
 23H  12.73000000   -9.73500000  -10.07300000
 24H  13.57000000  -10.91400000   -9.08500000
 25H  11.82400000  -10.73000000   -8.95800000
 26C  10.63500000  -17.66700000   -9.63800000
 27H  10.42500000  -18.14800000  -10.59400000
 28C   9.32200000  -17.82400000   -8.85000000
 29O   8.44600000  -18.59600000   -9.23100000
 30C  11.85000000  -18.46100000   -9.07300000
 31H  12.67900000  -18.38800000   -9.77700000
 32H  11.58400000  -19.51900000   -9.05700000
 33C  12.34800000  -18.12800000   -7.66200000
 34O  11.63600000  -18.31400000   -6.68200000
 35N   9.15900000  -16.99400000   -7.82000000
 36H   9.97600000  -16.44000000   -7.57700000
 37N  13.57300000  -17.63100000   -7.53200000
 38H  14.15700000  -17.47000000   -8.33800000
 39H  13.90400000  -17.39300000   -6.61000000
 40C   7.91700000  -16.59300000   -7.18800000
 41H   7.17100000  -16.53900000   -7.98100000
 42C   8.12600000  -15.13200000   -6.75300000
 43O   8.58800000  -14.32500000   -7.56700000
 44C   7.45500000  -17.63700000   -6.13700000
 45H   8.13900000  -17.64100000   -5.28900000
 46H   7.47900000  -18.64500000   -6.55300000
 47S   5.76100000  -17.29100000   -5.53200000
 48N   7.75900000  -14.75500000   -5.52100000
 49H   7.33100000  -15.46600000   -4.94500000
 50C   7.84900000  -13.41900000   -4.90300000
 51H   7.63700000  -13.58400000   -3.84500000
 52C   4.79400000  -17.78100000   -6.98300000
 53O   5.31900000  -18.26600000   -7.98300000
 54C   6.75400000  -12.45200000   -5.39100000
 55O   6.36800000  -11.52400000   -4.68400000
 56C   9.26200000  -12.81100000   -4.98000000
 57H   9.52900000  -12.52000000   -5.99400000
 58H  10.01500000  -13.51100000   -4.61600000
 59H   9.32600000  -11.91300000   -4.36400000
 60N   6.23100000  -12.71300000   -6.58600000
 61H   6.63300000  -13.48600000   -7.10100000
 62C   5.06800000  -12.13300000   -7.22300000
 63H   4.39300000  -11.78600000   -6.43700000
 64C   4.34200000  -13.27600000   -7.94400000
 65O   4.98800000  -14.22300000   -8.40500000
 66C   5.45000000  -10.90600000   -8.08100000
 67H   5.68700000  -10.08200000   -7.40700000
 68H   4.57100000  -10.57200000   -8.63400000
 69C   6.59600000  -11.05600000   -9.08100000
 70N   3.00400000  -13.22600000   -7.97800000
 71H   2.56500000  -12.44000000   -7.52100000
 72C   6.46800000  -11.89100000  -10.21300000
 73H   5.57000000  -12.47000000  -10.37200000
 74C   7.77500000  -10.29500000   -8.91600000
 75H   7.88200000   -9.63800000   -8.06300000
 76C   2.08200000  -14.21600000   -8.54300000
 77H   1.10200000  -13.74300000   -8.47200000
 78C   7.49800000  -11.95200000  -11.17200000
 79H   7.38000000  -12.57800000  -12.04600000
 80C   8.79400000  -10.33800000   -9.88700000
 81H   9.66700000   -9.70800000   -9.77900000
 82C   2.02900000  -15.45100000   -7.63000000
 83O   1.11000000  -15.59700000   -6.82700000
 84C   2.36100000  -14.49400000  -10.05000000
 85H   3.31200000  -15.01000000  -10.17600000
 86H   2.50200000  -13.52900000  -10.53900000
 87C   8.65300000  -11.16400000  -11.01900000
 88H   9.41100000  -11.17300000  -11.78900000
 89C   1.27900000  -15.25700000  -10.85600000
 90H   1.58500000  -15.16900000  -11.90000000
 91N   3.04600000  -16.31000000   -7.73200000
 92H   3.78900000  -16.04200000   -8.36900000
 93C   1.22600000  -16.76700000  -10.57000000
 94H   0.61600000  -17.28000000  -11.31500000
 95H   0.78900000  -16.98800000   -9.59700000
 96H   2.22000000  -17.21500000  -10.60300000
 97C  -0.11500000  -14.61600000  -10.76300000
 98H  -0.55300000  -14.73700000   -9.77100000
 99H  -0.80200000  -15.07400000  -11.47600000
100H  -0.07800000  -13.55000000  -10.98700000
101C   3.29500000  -17.48200000   -6.90300000
102H   3.05400000  -17.20100000   -5.87700000
103C   2.37100000  -18.65600000   -7.32500000
104H   2.62500000  -18.97100000   -8.33900000
105H   1.35100000  -18.27900000   -7.39600000
106C   2.32400000  -19.89900000   -6.40000000
107H   1.49800000  -20.50500000   -6.77600000
108C   3.56800000  -20.80000000   -6.48700000
109H   3.83500000  -21.01000000   -7.52400000
110H   3.38400000  -21.75900000   -6.00000000
111H   4.43500000  -20.35800000   -5.99900000
112C   1.97700000  -19.55500000   -4.94300000
113H   1.08500000  -18.93000000   -4.88600000
114H   2.79000000  -19.02600000   -4.44400000
115H   1.77800000  -20.45900000   -4.36600000

This molecule has a lot of rotatable bonds and some rings. In principle, you can manually provide these flexibility information for geom, but you can use to topgen instead:

$ topgen aip.xyz

There are 3 files related to this work:

  • aip.gjf: The molecule in GJF format containing bonding information for reference.

  • aip-bonding.xyz: The molecule in XYZ format containing bonding information. You can find the following after coordinates:

aip-bonding.xyz
 1H   2.79000000  -19.02600000   -4.44400000
 2H   1.77800000  -20.45900000   -4.36600000
 3
 41 2 1.0 3 1.0 4 1.0 5 1.0
 52
 63
 74
 85 6 1.0 7 1.0 9 1.0
 96
107 8 2.0 18 1.0
118
129 10 1.0 11 1.0 14 1.0
1310
14.. omitted ..
15110 111 1.0 112 1.0 113 1.0
16111
17112
18113

This bonding information is completely compatible with GaussView. Since there are 113 atoms, there will be 113 lines of bonding information for every atom. 1 2 1.0 3 1.0 4 1.0 5 1.0 means that atom 1 forms a bond with order 1.0 with atom 2, 3, 4, and 5. Note that only the bonding atoms with indices larger than this atom are written. For example, 5 6 1.0 7 1.0 9 1.0 `` indicates that atom ``5 bonds with atom 6, 7, and 9. Its bonding with atom 1 is not written here but in the line 1 2 1.0 3 1.0 4 1.0 5 1.0. Only the bonds with order 1.0 are treated as flexible bonds. The ones with order 2.0, 3.0, etc are treated as non-rotatable.

You must provide an XYZ file with bonding information to geom if you want ABCluster to explicitly consider flexibility.

  • aip-cycles.txt: It lists cycles in this molecule.

aip-cycles.txt
138 40 46 48 52 58 60 62 68 74 80 89 99 50 45 42
267 70 76 85 78 72

This means that there are two rings in this molecule: a 16-membered ring composed of atom 38, 40, …, 45, and 42 and a 6-membered ring composed of atom 67, …, 72.

These lines can be directly used in the input of geom.

Tip

Note that topgen does not detect if a ring is rigid or not. For the above example, ring 67 70 76 85 78 72 is a benzene so it should not be deformed during conformation search.

7.1.2. Set up Flexibility in Input Files

You need to set flexibility in component block of an input file:

geom.inp
1components
2  aip-bonding.xyz 1
3  fix 0 0 0 0 0 0
4  ****
5  water.xyz 100
6  random -5 -5 -5 5 5 5
7  ****
8end

As described in Control the Cluster Shape, this means that we fix molecule AIP at position (0, 0, 0) without rotation, and randomly put 100 water in the box defined by (-5, -5, -5) and (5, 5, 5). Note that now AIP is still rigid since you do not assign flexibility.

You can set flexibility between the placement way and ****:

geom.inp
 1components
 2  aip-bonding.xyz 1
 3  fix 0 0 0 0 0 0
 4  nH nH
 5  38 40 46 48 52 58 60 62 68 74 80 89 99 50 45 42
 6  ****
 7  water.xyz 100
 8  random -5 -5 -5 5 5 5
 9  ****
10end

Here nH nH and 38 40 46 48 52 58 60 62 68 74 80 89 99 50 45 42 are the flexible degrees of freedom to be sampled. There can be arbitrary lines of flexibility:

  • 5 or more atomic indices: defines a flexible ring. This ring will be deformed during global optimization by geom. Say, 38 40 46 48 52 58 60 62 68 74 80 89 99 50 45 42 is a ring defined by atom 38, …, 42.

  • 2 characters: defines flexible bonds. These bonds will be rotated during global optimization by geom. Note that there are two useful wildcards: * (any atom) and nH (any non-hydrogen atom). Now bonds can be defined in many forms:

    • 23 49 The bond between atom 23 and 49.

    • C C All bonds between C and C.

    • C N All bonds between C and N.

    • nH nH All bonds between non-hydrogen atoms.

    • * * All bonds.

    • 15 * All bonds between atom 15 and any atoms.

    • 15 nH All bonds between atom 15 and any non-hydrogen atoms.

7.1.3. Some Examples

geom.inp
1components
2  aip-bonding.xyz 1
3  fix 0 0 0 0 0 0
4  nH nH
5  ****
6end

All bonds between non-hydrogen atoms are rotatable.

geom.inp
1components
2  aip-bonding.xyz 1
3  fix 0 0 0 0 0 0
4  18 19
5  45 50
6  64 nH
7  67 70 76 85 78 72
8  ****
9end

Bonds between 18 and 19, 45 and 50, 64 and any non-hydrogen atoms, are rotatable. The ring 67 70 76 85 78 72 is deformable.

geom.inp
1components
2  aip-bonding.xyz 1
3  fix 0 0 0 0 0 0
4  C C
5  ****
6end

All bonds between C and C are rotatable.

geom.inp
1components
2  aip-bonding.xyz 1
3  fix 0 0 0 0 0 0
4  nH nH
5  38 40 46 48 52 58 60 62 68 74 80 89 99 50 45 42
6  67 70 76 85 78 72
7  ****
8end

All bonds between non-hydrogen atoms are rotatable. The ring 67 70 76 85 78 72 and 38 40 46 48 52 58 60 62 68 74 80 89 99 50 45 42 are deformable.

Now you can use geom to do conformation search and global optimization of clusters of highly flexible units!

Tip

As mentioned above, ring 67 70 76 85 78 72 is a benzene so it should NOT be deformed during conformation search. But if you write here, it will be forcefully deformed, although during optimization it may probably go back to a plane.