求助有关gromacs求扩散系数

少年爱吃地三鲜

请问老师， GROMACS 内部求扩散系数是 通过纵坐标之差与横坐标之差的比来表示的吗？（除以6）

liuyuje714

就是用10%~90%这个时间区域内的数据做最小二乘拟合得到斜率a，然后根据体系扩散类型，比如三维自由扩散，就除以2*DIM=6，然后乘1000转换单位，下面给你一段awk代码做拟合直线

1. #!/bin/bash
   
2. 
   
3. DATA=$1
   
4. awk -v data=$DATA '
   
5.     BEGIN {
   
6.         sumxy = 0; sumxx = 0; sumx = 0; sumy = 0; n = 0;
   
7.         ss_t = 0; ss_r = 0;
   
8.         while(getline<data) {
   
9.             if($0!~/#|@/)
   
10.             {
    
11.                 sumxy += $1 * $2
    
12.                 sumxx += $1 * $1
    
13.                 sumx  += $1
    
14.                 sumy  += $2
    
15.                 n++
    
16.             }
    
17.         }
    
18.         close(data)
    
19.         k = (sumxy - n * (sumx/n) * (sumy/n))/(sumxx - n * (sumx/n)^2)
    
20.         b = (sumy/n) - (k * (sumx/n))
    
21.         printf("Evaluation of equation: y = (%.14f) * x + (%.14f)\n", k, b)
    
22.     }
    
23.    
    
24.     {
    
25.     ss_t += ($2 - sumy/n)^2
    
26.     ss_r += ($2 - (k * $1 + b))^2
    
27.     }
    
28.    
    
29.     END {
    
30.         r2 = 1 - (ss_r / ss_t)
    
31.         printf("\t\t\tR^2 = %.6f\n", r2)
    
32.     }
    
33. ' $DATA
    

复制代码

少年爱吃地三鲜

liuyuje714 发表于 2020-10-28 16:55
就是用10%~90%这个时间区域内的数据做最小二乘拟合得到斜率a，然后根据体系扩散类型，比如三维自由扩散，就 ...

感谢您， 我看一下

少年爱吃地三鲜

liuyuje714 发表于 2020-10-28 16:55
就是用10%~90%这个时间区域内的数据做最小二乘拟合得到斜率a，然后根据体系扩散类型，比如三维自由扩散，就 ...

您好，我想请教假如我想求扩散系数随着时间的变化情况，这个D 是否可以安装公式，先求某一点的斜率，然后除以6，这个点的斜率 我用相邻两个点的坐标求近似，请问这样是否是可行的？

liuyuje714

少年爱吃地三鲜 发表于 2020-10-28 21:55
您好，我想请教假如我想求扩散系数随着时间的变化情况，这个D 是否可以安装公式，先求某一点的斜率，然后 ...

这就涉及到一个拟合区域的选择，一般开头和结尾的数据都不应该选。即计算动扩散系数，就是Msd曲线对关联时间的导数dMSD(t)/6dt，即你所说斜率，一般选择成水平直线的部分。

少年爱吃地三鲜

liuyuje714 发表于 2020-10-28 22:25
这就涉及到一个拟合区域的选择，一般开头和结尾的数据都不应该选。即计算动扩散系数，就是Msd曲线对关联 ...

我明白了！感谢您!

ocbrother

liuyuje714 发表于 2020-10-28 22:25
这就涉及到一个拟合区域的选择，一般开头和结尾的数据都不应该选。即计算动扩散系数，就是Msd曲线对关联 ...

你好，请问计算扩散系数选水平直线的部分的话，斜率不应该趋近于0了吗。。

sobereva

ocbrother 发表于 2020-12-14 21:40
你好，请问计算扩散系数选水平直线的部分的话，斜率不应该趋近于0了吗。。

不是水平直线，是直线

ocbrother

sobereva 发表于 2020-12-16 04:51
不是水平直线，是直线

请问sob老师，MSD-t图上有很多斜率不同的直线段时，该取哪一段作为扩散系数呢？

liuyuje714

ocbrother 发表于 2020-12-16 09:00
请问sob老师，MSD-t图上有很多斜率不同的直线段时，该取哪一段作为扩散系数呢？

我说的是msd/6t的水平直线，我也没说是0啊，你自己想啊，msd如果是成线性那他的斜率msd/t必然是水平直线啊，扩散系数就是斜率/6。我没用过MS，所以不知道具体情况。

ocbrother

liuyuje714 发表于 2020-12-16 09:08
我说的是msd/6t的水平直线，我也没说是0啊，你自己想啊，msd如果是成线性那他的斜率msd/t必然是水平直线 ...

对不起 是我误会你的意思了。不好意思

sobereva

ocbrother 发表于 2020-12-16 09:00
请问sob老师，MSD-t图上有很多斜率不同的直线段时，该取哪一段作为扩散系数呢？

没具体图不好说。一般取相对靠中间的，最长的一段

ocbrother

sobereva 发表于 2020-12-17 09:40
没具体图不好说。一般取相对靠中间的，最长的一段

好的 谢谢老师