marrink-lab / polyply_1.0

Generate input parameters and coordinates for atomistic and coarse-grained simulations of polymers, ssDNA, and carbohydrates
Apache License 2.0
122 stars 21 forks source link

Add PTMA AA FF to polyply library [10.1021/acs.macromol.3c00141] #276

Closed ricalessandri closed 1 year ago

ricalessandri commented 1 year ago

PR for contributing the OPLS-based AA FF developed for PTMA in https://arxiv.org/abs/2209.02072.

I will shortly add a test. Also open a separate PR to contribute CG model(s).

fgrunewald commented 1 year ago

@ricalessandri if you register here you can get a free T-shirt if you make 4 PRs in Oktober for the Hacktoberfest.

ricalessandri commented 1 year ago

I have not looked into that (sorry, I remember about our discussions but did not have the bandwidth for that). Maybe this is a good chance. I will try to produce a plot asap.

I am however unsure whether I would really want to introduce things in this model at this point since this is the model that I have been using for the simulations.

fgrunewald commented 1 year ago

Sure thing we just have to know if it gives an atactic or isotactic polymer. Either way is fine for me.

ricalessandri commented 1 year ago

Finally getting to close this, some things to fix before merging so I don't forget:

ricalessandri commented 1 year ago

Alright, I'm finally at it. This is one of the 3 things left to close this PR: tacticity. @fgrunewald we're getting atactic PTMAs, see:

Screen Shot 2023-08-21 at 8 50 45 AM
Average enantiomer fractions over 100 chains of 30-mers:
enantiomer #1 : 0.39399930000000005
enantiomer #2 : 0.6059994800000001

With a 40-60% ratio. I guess this would converge to 50-50 with more (and/or longer) chains. Also, chains are not syndiotactic. I checked this by printing the dihedrals. See for yourself for the first 10 chains:

  0 36.5249         |  0 36.0532    |  0 34.7708    |  0 -36.5657   |  0 35.2323    |  0 36.9282    |  0 36.1829    |  0 36.1268    |  0 -36.0353
  1 34.0585         |  1 32.8891    |  1 35.7504    |  1 34.8827    |  1 38.3678    |  1 -41.1577   |  1 40.1183    |  1 -35.7562   |  1 37.7092
  2 -37.2684        |  2 34.5314    |  2 -36.087    |  2 36.5295    |  2 38.1303    |  2 32.8347    |  2 37.3767    |  2 -35.1527   |  2 -39.0089
  3 36.8553         |  3 34.2684    |  3 -36.7091   |  3 36.5085    |  3 -36.8466   |  3 -40.7319   |  3 35.2229    |  3 40.2213    |  3 40.0856
  4 -33.6193        |  4 -32.8512   |  4 -42.2365   |  4 -37.6941   |  4 -34.7345   |  4 -37.4294   |  4 -36.0298   |  4 -36.4733   |  4 35.8881
  5 37.2838         |  5 33.4698    |  5 -32.5688   |  5 37.3802    |  5 35.678     |  5 36.2165    |  5 37.1842    |  5 37.5854    |  5 -32.3946
  6 35.3496         |  6 -36.0319   |  6 38.4035    |  6 39.7357    |  6 -41.5397   |  6 38.9863    |  6 40.6138    |  6 36.6165    |  6 37.2975
  7 32.874          |  7 -36.0374   |  7 36.9036    |  7 36.3968    |  7 38.5483    |  7 -39.3083   |  7 -35.7891   |  7 37.9968    |  7 34.7671
  8 43.9699         |  8 35.4649    |  8 -37.1433   |  8 -35.7173   |  8 37.0971    |  8 36.0991    |  8 -33.8231   |  8 31.2062    |  8 36.5176
  9 35.8624         |  9 37.4412    |  9 38.3855    |  9 -39.6864   |  9 -41.669    |  9 -38.5243   |  9 38.1871    |  9 40.1058    |  9 -37.2323
  10 37.6172        |  10 39.4304   |  10 34.4706   |  10 -39.4187  |  10 35.063    |  10 39.1183   |  10 33.5543   |  10 -34.6115  |  10 38.5975
  11 41.7453        |  11 33.5612   |  11 35.5474   |  11 34.1137   |  11 -36.8871  |  11 35.875    |  11 -38.4814  |  11 37.2663   |  11 41.8513
  12 39.0329        |  12 -34.4066  |  12 36.363    |  12 31.7394   |  12 36.6417   |  12 -35.4788  |  12 -33.0636  |  12 -35.5622  |  12 -36.0089
  13 35.9318        |  13 39.8196   |  13 -39.2512  |  13 35.2016   |  13 -33.0821  |  13 -36.4384  |  13 39.1468   |  13 37.0456   |  13 36.0463
  14 -34.3356       |  14 -36.6234  |  14 37.2354   |  14 39.4331   |  14 41.5223   |  14 37.71     |  14 -38.5657  |  14 37.4491   |  14 37.891
  15 -34.7859       |  15 -33.715   |  15 -36.5078  |  15 -38.4706  |  15 31.7546   |  15 37.0417   |  15 36.2768   |  15 -33.8178  |  15 -37.1269
  16 33.7626        |  16 -36.9238  |  16 33.5581   |  16 40.7569   |  16 36.9485   |  16 34.397    |  16 35.0809   |  16 36.8588   |  16 37.7889
  17 39.1062        |  17 38.5846   |  17 33.7556   |  17 38.6473   |  17 37.9576   |  17 36.227    |  17 33.3178   |  17 38.425    |  17 36.5993
  18 -36.0831       |  18 33.9731   |  18 -39.3478  |  18 -36.6049  |  18 -35.2307  |  18 -35.9973  |  18 33.6778   |  18 38.6039   |  18 34.586
  19 -32.5932       |  19 -34.7784  |  19 -34.793   |  19 -34.527   |  19 38.7232   |  19 -33.607   |  19 33.8859   |  19 35.1522   |  19 -35.5263
  20 -39.7191       |  20 38.3118   |  20 35.6131   |  20 36.365    |  20 37.1533   |  20 39.9286   |  20 37.0797   |  20 -36.2315  |  20 36.6258
  21 39.7267        |  21 39.9699   |  21 38.0917   |  21 37.1777   |  21 -42.0094  |  21 36.3877   |  21 -36.6893  |  21 -36.845   |  21 -34.415
  22 34.2142        |  22 36.3733   |  22 35.9077   |  22 -35.3277  |  22 36.8763   |  22 36.5706   |  22 36.445    |  22 -33.9113  |  22 40.454
  23 -39.275        |  23 37.1109   |  23 -37.8667  |  23 38.5636   |  23 -35.4948  |  23 -39.1039  |  23 36.604    |  23 34.715    |  23 -35.9212
  24 36.3931        |  24 37.9223   |  24 37.0983   |  24 37.749    |  24 -33.4265  |  24 35.8632   |  24 -37.5184  |  24 36.9496   |  24 -38.328
  25 -34.1988       |  25 -39.3313  |  25 38.5158   |  25 37.0989   |  25 -37.5111  |  25 37.9106   |  25 36.1155   |  25 -38.0322  |  25 38.8457
  26 34.4993        |  26 29.5963   |  26 -35.9228  |  26 37.5261   |  26 36.3587   |  26 37.3568   |  26 -37.3268  |  26 -35.1047  |  26 39.0637
  27 33.1792        |  27 35.752    |  27 -37.0754  |  27 -38.6848  |  27 37.4816   |  27 -35.7833  |  27 38.576    |  27 -34.3623  |  27 37.3753
  28 -34.6558       |  28 42.5101   |  28 -40.0315  |  28 37.6345   |  28 -38.7036  |  28 35.2331   |  28 -39.6963  |  28 37.8905   |  28 32.9077
  29 31.8568        |  29 31.8861   |  29 40.188    |  29 -37.5955  |  29 33.4013   |  29 -38.7681  |  29 -36.739   |  29 37.571    |  29 38.1437

Not sure if you want this more quantitative. In case, what's a smart metric for this?

ricalessandri commented 1 year ago

I resolved the previous issues. I also added a test: 2e207a3. It seems to be working.

I guess this is ready to be merged, @fgrunewald .

ricalessandri commented 1 year ago

One test fails but I think it's not related to this .ff file. EDIT: I re-triggered and all is fine.

fgrunewald commented 1 year ago

@ricalessandri you can also consider adding a line on the news section of the readme file to promote your paper

ricalessandri commented 1 year ago

Thanks for the nice ideas!

I've added the warning, it works: https://github.com/marrink-lab/polyply_1.0/blob/006e37ed44d2f8730e806cd66d5de8b71c4859f0/polyply/data/oplsaaLigParGen/links_PTMA.ff#L164-L174

I didn't include that [ non-edges ] sections of P3HT, which would have looked like this:

[ link ]
[ atoms ]
C01 {"resname": "PTMA"}
[ warning ]
You should patch PTMA with a CH3-terminal residues. Use -seq CH3a:1 PTMA:4 CH3b:1
[ non-edges ]
-C1 C01
-C17 C01

[ link ]
[ atoms ]
C17 {"resname": "PTMA"}
[ warning ]
You should patch PTMA with a CH3-terminal residues. Use -seq CH3a:1 PTMA:4 CH3b:1
[ non-edges ]
C17 +C1
C17 +C01

because it didn't make any difference in the resulting itp and I wasn't sure about their purpose. Let me know if it should be there.

As a side note, I played with the behavior and the warnings will be printed twice for each residue (i.e., 8 times for PTMA:4). While the twice is fine because it's one for CH3a and one for CH3b, printing it for each residue seems unnecessary. Also, if I use: PTMA:4 CH3b:1 or CH3a:1 PTMA:4 I was expecting polyply to print only one warning, but it still prints two. Let me know if there's something wrong with my definitions. Otherwise, these are anyway minor things because we want the user to use CH3a:1 PTMA:4 CH3b:1 anyway.

I will include the news piece in the other PR later today.

fgrunewald commented 1 year ago

@ricalessandri the non-edges are specifically to reduce the number of warnings being printed and in principle, you could even use them to specifically target only one end, like if you miss CH3b it would only do 1 warning. However, there is a small syntax mistake in the non-edges of the first residue. It should first mention the atom that is also in the atoms directive. So it should be:

[ link ]
[ atoms ]
C01 {"resname": "PTMA"}
[ warning ]
You should patch PTMA with a CH3-terminal residues. Use -seq CH3a:1 PTMA:4 CH3b:1
[ non-edges ]
C01 -C1
C01 -C17

If you adjust this, the PR should be ready to go.

ricalessandri commented 1 year ago

@fgrunewald done! I implemented and tested and works as expected (1 warning for CH3a if I miss that, 1 warning for CH3b if I miss that, of 2 warnings if I miss both) - fantastic!