Bug in CCFD stress profile calculation of GROMACS-LS?

In Martini simulations (15 us, last 5 us analyzed) of an asymmetric GM1+POPC lipid bilayer (5+95 upper, 0+100 lower leaflet) strong stress fluctuations arising from the dihedral contributions of CCF decomposition are found. It is possible that these are due to a slight bug in the GROMACS-LS (http://mdstress.org) code, because: 1) The fluctuations do not show converge, when more data is added, see the attached plot ‘effect_of_sampling_rate.pdf’. 2) The fluctuations (but no other features of the stress profiles) visibly depend on the chosen frame of reference, see ‘effect_of_centering_schemes.pdf’. 3) The stress profiles calculated along the coordinates that lay in the membrane plane (x and y), are not constant in CCFD, see ‘profiles_along_x_and_y.pdf’. To assist debugging the code, this repository contains: <strong>5to95_0to100_0.trr</strong> --- 1st 33.3% of a 5 us trajectory with snapshots every 100 ps<br> <strong>5to95_0to100_0.trr</strong> --- 2nd 33.3% of a 5 us trajectory with snapshots every 100 ps<br> <strong>5to95_0to100_0.trr</strong> --- 3rd 33.3% of a 5 us trajectory with snapshots every 100 ps <strong>5to95_0to100.tpr</strong> --- tpr used for production (Gromacs 5.1.1)<br> <strong>5to95_0to100_rerun.tpr</strong> --- tpr used for stress analysis and centering <strong>5to95_0to100.gro</strong> --- gro file after 15 us (10 us relaxation, 5 us production) <strong>index.ndx</strong> --- index file used in production<br> <strong>5to95_0to100.ndx</strong> --- index file used in analysis <strong>system.top</strong> --- topology file<br> <strong>martini_v2.0_ions.itp</strong> --- martini FF file (system.top expects this to be in folder toppar/)<br> <strong>martini_v2.0_lipids_all_201506.itp</strong> --- martini FF file (system.top expects this to be in folder toppar/)<br> <strong>martini_v2.2.itp</strong> --- martini FF file (system.top expects this to be in folder toppar/) <strong>martini_straight_GM1_saveVels.mdp</strong> --- run input file used for production (Gromacs 5.1.1)<br> <strong>martini_straight_GM1_rerunForP.mdp</strong> --- run input file used for analysis and centering <strong>effect_of_sampling_rate.pdf</strong> --- plot of results showing that fluctuations in CCFD do not appear to converge when data is added<br> <strong>profiles_along_x_and_y.pdf</strong> --- plot of results showing that in CCFD pressure profiles along the membrane directions are not constant<br> <strong>effect_of_centering_schemes.pdf</strong> --- plot of results showing that choice of center of mass (here three possibilities are shown: CoM of POPC GL1 beads, CoM of POPCs, and CoM of all lipids) visibly affects CCFD fluctuations, but not other features of the stress profile<br> <strong>effect_of_centering_schemes_2.pdf</strong> --- plot of results (similar to 'effect_of_centering_schemes.pdf', but the center of the bilayer is at 3.7 nm instead at 0.0 / 7.4 nm)<br> <strong>individual_components.pdf</strong> --- plot of results showing that the fluctuations in CCFD arise from dihedrals<br> <strong>CCFD_gridsize.pdf</strong> --- plot of results showing that decreasing gridsize has no clear effect on the CCFD fluctuations本数据集针对非对称GM1+POPC脂质双层（上层叶瓣含5% GM1、95% POPC，下层叶瓣含0% GM1、100% POPC）开展了Martini模拟（总时长15 μs，仅分析最后5 μs的模拟数据），观测到由CCF分解（CCF decomposition，下文简称CCFD）的二面角贡献项引发的显著应力涨落现象。 推测该现象可能源于GROMACS-LS（http://mdstress.org）代码存在细微漏洞，依据如下： 1. 随着新增更多模拟数据，应力涨落并未呈现收敛趋势，详见附件绘图`effect_of_sampling_rate.pdf`（采样率影响.pdf）； 2. 应力涨落（而非应力剖面的其他特征）明显依赖于所选参考帧，详见`effect_of_centering_schemes.pdf`（质心校正方案影响.pdf）； 3. 沿膜平面坐标（x与y方向）计算得到的应力剖面，其CCFD值并非恒定，详见`profiles_along_x_and_y.pdf`（x/y方向剖面绘图.pdf）。 为协助代码调试，本仓库包含以下文件： - `5to95_0to100_0.trr`：5 μs模拟轨迹的前33.3%片段，每100 ps保存一个快照 - `5to95_0to100_0.trr`：5 μs模拟轨迹的中间33.3%片段，每100 ps保存一个快照 - `5to95_0to100_0.trr`：5 μs模拟轨迹的后33.3%片段，每100 ps保存一个快照 - `5to95_0to100.tpr`：生产模拟所用的TPR文件（适配GROMACS 5.1.1版本） - `5to95_0to100_rerun.tpr`：用于应力分析与坐标校正的TPR文件 - `5to95_0to100.gro`：15 μs模拟结束后的GRO结构文件（其中前10 μs为弛豫阶段，后5 μs为生产模拟阶段） - `index.ndx`：生产模拟所用的索引文件 - `5to95_0to100.ndx`：分析阶段所用的索引文件 - `system.top`：系统拓扑文件 - `martini_v2.0_ions.itp`：Martini力场离子参数文件（系统拓扑文件`system.top`要求该文件存放于`toppar/`目录下） - `martini_v2.0_lipids_all_201506.itp`：Martini力场脂质参数文件（系统拓扑文件要求该文件存放于`toppar/`目录下） - `martini_v2.2.itp`：Martini力场通用参数文件（系统拓扑文件要求该文件存放于`toppar/`目录下） - `martini_straight_GM1_saveVels.mdp`：生产模拟所用的运行输入文件（适配GROMACS 5.1.1版本） - `martini_straight_GM1_rerunForP.mdp`：用于应力分析与坐标校正的运行输入文件 - `effect_of_sampling_rate.pdf`：展示随模拟数据量增加，CCFD涨落未呈现收敛趋势的结果绘图 - `profiles_along_x_and_y.pdf`：展示沿膜平面方向的CCFD压力剖面并非恒定的结果绘图 - `effect_of_centering_schemes.pdf`：展示质心（center of mass, CoM）选择对CCFD涨落存在显著影响，但不影响应力剖面其他特征的结果绘图（此处展示三种质心方案：POPC GL1珠子的质心、所有POPC的质心、所有脂质的质心） - `effect_of_centering_schemes_2.pdf`：与`effect_of_centering_schemes.pdf`类似的结果绘图，但双层膜的中心位置为3.7 nm，而非0.0 nm或7.4 nm - `individual_components.pdf`：展示CCFD涨落源自二面角贡献项的结果绘图 - `CCFD_gridsize.pdf`：展示网格尺寸减小对CCFD涨落无显著影响的结果绘图