|
|
各位老师好，我是模拟离子穿越细胞膜的过程，膜上下离子浓度差异的扩散跨膜运输。但是由于周期性边界条件XYZ的设置，出现了膜上的部分离子穿过周期性镜像从最上面穿过，由最下面出现，从下往上运动了，导致了不合适的轨迹出现。因此我想在Z方向设定wall，pbc=xy来防止这样的现象出现，但是作为新手不是很会设置这个，导致体系出现了崩溃，所以想请教各位老师wall这里如何设置？并有以下疑问：
1.我尝试过在体系上下设置足够大的真空来防止这个现象出现，那这样是不是就只能使用NVT进行预平衡和成品模拟了？（我在纠结根据一般规则好像都是要NVT预平衡+NPT预平衡+成品，留有真空部分的话貌似NPT就不能用了，但是这样是不是就不能对体系中的膜蛋白充分预平衡了，我只是想观察和分析离子从上而下跨膜扩散运输的过程，以及和蛋白的结合位点，不用常规的预平衡中的NPT平衡膜蛋白体系的话，担心会有审稿人的质疑，不知道是不是自己跨入了误区，恳请各位老师给予帮助）
2.请问wall具体如何设置较为合适，我是在成品MD模拟中设置的wall，NVT预平衡+NPT预平衡阶段没有设置wall，只是把离子都设置了位置限制，成品MD阶段解除了位置限制，以下是我的MD.mdp文件设置：
;{ 预处理, 使用 C++ 语法
;===============================================================================
include = ; 含引用文件的目录, 拓扑文件中可引用其中的文件, 可多项
; 例: -I/home/joe/joy -I/home/tom/tim -IC:/GMX/top
define = -DFLEXIBLE ; 预定义, 默认无, 可多项, 区分大小写
; -DPOSRES: 使用位置限制文件进行位置限制动力学模拟
; -DFLEXIBLE: 启用柔性水, steep效果更好, cg, l-bfgs或简正分析须开启
;}==============================================================================
;{ MD运行控制
;===============================================================================
integrator = md ; 积分方法, md:蛙跳; sd:随机; steep/cg/lbfgs:能量最小化
dt = 0.002 ; 积分步长(ps), EM不用
nsteps = 250000 ; 最大积分步数, 默认0, -1:无限制
comm-mode = Linear ; 移除质心运动的方式, None:无; Linear:平动; Angular:平动转动;
; Linear-acceleration-correction:线性加速度校正, 牵引情况
nstcomm = 100 ; 移除质心运行的频率(步)
comm-grps = SOLU_MEMB SOLV ; 移除质心运动的组, 可多个, 默认整个体系
tinit = 0 ; 起始时间(ps), EM不用
init-step = 0 ; 起始步数, 对非平衡模拟, 精确重启或重做某部分模拟时, 设定为重启步编号
simulation-part = 1 ; 检查点时自动更新的部分编号(保持文件分开)
;}==============================================================================
;{ 输出控制
;===============================================================================
nstxout = 1000 ; trr坐标的输出频率(步)
nstvout = 1000 ; 速度
nstfout = 1000 ; 力
nstxout-compressed = 1000 ; xtc压缩坐标的输出频率
nstlog = 1000 ; 日志文件输出频率
nstenergy = 1000 ; 能量文件输出频率
nstcalcenergy = 1000 ; 计算能量的频率, 最好为 nstlist 倍数
energygrps = ; 输出到能量文件的组, 可使用多个, 默认所有
compressed-x-grps = ; 输出xtc压缩坐标的组, 可使用多个, 默认所有
compressed-x-precision = 1000 ; xtc坐标的精度(1000表示1/1000, 三位小数)
;}==============================================================================
;{ 邻区搜索
;===============================================================================
cutoff-scheme = Verlet ; 截断方式, Verlet:粒子截断; Group:电荷组
ns-type = Grid ; 邻区搜索算法, Grid:较快; Simple:仅与Group联用
nstlist = 20 ; 邻区列表更新频率, 0:真空模拟; -1:自动
rlist = 1.2 ; 邻区列表截断距离(nm)
rvdw = 1.2 ; 范德华截断半径
rcoulomb = 1.2 ; 静电截断半径, 不超过最小盒子边长一半
nstcalclr = -1 ; 长程邻区列表的计算频率
rlistlong = -1 ; 切换势能函数的长程邻区列表截断距离(nm)
verlet-buffer-tolerance = 0.005 ; Verlet缓冲的能量误差(kJ/mol-ps-atom), -1:使用rlist
;}==============================================================================
;{ 周期性
;===============================================================================
pbc = XY ; 周期性边界条件, XYZ; XY; No:忽略盒子, 截断与nstlist置零
periodic-molecules = No ; 周期性分子: No; Yes
;}==============================================================================
;{ 静电与范德华
;===============================================================================
vdwtype = Cut-off ; 范德华计算方法
coulombtype = PME ; 静电计算方法
DispCorr = No ; 长程色散校正, No:无; Ener:能量; EnerPres:能量和压力
rvdw-switch = 0 ; 切换距离
rcoulomb-switch = 0
vdw-modifier = Force-switch ; 修正方法
epsilon-r = 1 ; 介质的相对介电常数, 0:无穷大
epsilon-rf = 0 ; 反应场的相对介电常数, 0:无穷大
;}==============================================================================
;{ 初始速度
;===============================================================================
gen-vel = No ; 产生方式, Yes:随机; No:使用gro文件中的值
gen-temp = 303.15 ; 随机速度对应的温度
gen-seed = -1 ; 随机数种子; -1:自动确定
;}==============================================================================
;{ 温度耦合
;===============================================================================
tcoupl = Nose-Hoover ; 耦合方法, No:无; V-rescale:快速; Nose-Hoover:精确
tc-grps = SOLU MEMB SOLV ; 温度耦合组, 可多个
tau-t = 1.0 1.0 1.0 ; 时间常数(ps)
ref-t = 303.15 303.15 303.15 ; 参考温度(K)
nsttcouple = -1 ; 耦合频率, -1:同nstlist
nh-chain-length = 10
print-nose-hoover-chain-variables = No
;}==============================================================================
;{ 压力耦合
;===============================================================================
pcoupl = Parrinello-Rahman ; 耦合方法, No:无, 盒子大小不变; Berendsen:快速; Parrinello-Rahman:精确
pcoupltype = semiIsotropic ; 耦合类型, Isotropic:各向同性;
; semiIsotropic:x/y方向各向同性, 与z方向不同, 膜模拟
; anIsotropic:各向异性, 盒子可能剧烈变形
; surface-tension:表面张力
tau-p = 5 ; 时间常数(ps)
compressibility = 4.5E-5 4.5E-5 ; 压缩率(1/bar)
ref-p = 1 1 ; 参考压力(bar)
nstpcouple = -1 ; 耦合频率, -1:同nstlist
refcoord-scaling = COM ; 缩放参考坐标: No:无; All:所有粒子; COM:质心
;}==============================================================================
;{ 键约束
;===============================================================================
constraints = h-bonds ; None:仅拓扑中指定的约束
; h-bonds:含氢键; h-angles:含氢键和键角
; all-bonds:所有键; all-angles:所有键和键角
constraint-algorithm = Lincs ; 约束算法, Lincs:不支持键角; Shake:慢且不稳定, 不支持EM和dd
continuation = Yes ; 是否约束初始构型, 并重置壳层
Shake-SOR = No ; 使用连续超弛豫方法减少shake迭代数目
shake-tol = 1E-4 ; shake的相对容差
lincs-order = 4 ; 约束耦合矩阵展开的最高阶数
lincs-iter = 1 ; lincs最终步的迭代数, 1:常规模拟; 2:NVE能量守恒; 4或8:含有约束的能量最小化
lincs-warnangle = 30 ; 如果某一步中键的旋转角度超过此值, 显示lincs警告
morse = No ; 将简谐键转换为Morse势
;}==============================================================================
;{ 墙
;===============================================================================
nwall = 2 ; 墙的数目: 1:z=0; 2:z=0, z=z-box. 只能使用pbc=xy
wall-atomtype = CAD CLA ; 墙的原子类型名称
wall-type = 9-3 ; LJ势类型: 9-3:体积分; 10-4:面积分; 12-6:直接
table = ; 用户自定义的势能, 使用与墙的z距离作为索引
wall-r-linpot = 2 ; 与墙的距离在此值以下时, 势能线性连续. 正值有助于解决原子超出墙的问题(nm)
wall-density = 450 450 ; 每面墙的原子数密度, 用于9-3(nm^-3)和10-4(nm^-2)类型
wall-ewald-zfac = 3 ; 盒子z方向Ewald的缩放因子, 最小值2. Ewald只能与nwall=2, ewald-geometry=3dc联用
ewald-geometry = 3dc ; Ewald几何结构
;}==============================================================================
|
|
发表于 Post on 2023-10-27 11:15:38
|
只看该作者 Only view this author
|倒序浏览 Reverse view
|阅读模式 Reading model
收藏 Add to favorites