求助NaCl电解质体系总偶极矩模的计算

wangyueda

请问下各位老师，利用gmx对NaCl电解质体系进行模拟（溶剂为水），体系平衡后，如果我想计算一段时间内体系总偶极矩的模，即|<M>|（尖括号表示系综平均），应该如何利用tcl编程处理呢？
此处之所以没有用gmx dipoles来处理主要是因为，体系中存在大量的 [(Na)(Cl)(H2O)n]结构，这种结构总体净电荷为0，结构中NaCl视为分子，对体系总偶极矩的贡献不可忽略。显然，这种结构的数目在模拟过程中是不断变化的，也就是说此结构中NaCl是一个动态选区，因此不好用gmx dipoles来处理。我想着用tcl语言来处理，但是有下面几个问题:

(1) 总偶极矩系综平均的问题：根据总偶极矩的定义 M = ∑ (q r)，以三维坐标的形式表示，假设有模拟轨迹有1000帧结构，那么每一帧都可以算出来一个对应的 M，即第一帧 M1（x1,y1,z1）第二帧 M2（x2,y2,z2）...第一千帧 M1000（x1000,y1000,z1000）。那么总偶极矩系综平均<M>是将M1到M1000直接做数值平均吗 ？ 即<M>=（1/1000）（x1+x2+...x1000,y1+y2+...y1000,z1+z2+...z1000）


(2) 单位转换问题：偶极矩的单位Debye = Bohr * e，tcl脚本 $sel get x y z得到的坐标单位是A，那么此处需要将A转换为Bohr，即乘以1.89再计算偶极矩吗？

sobereva

1 是

2 是

wangyueda

sobereva 发表于 2023-1-30 15:01
1 是

2 是

好的 谢谢卢老师 ps 后续发现可以用measure dipole -debye来直接计算偶极矩，没必要写一堆循环了