CASEA的奇思妙想

APBS 输入文件

发表于 分类于
本文字数： 12k 阅读时长 ≈ 11 分钟

APBS 输入文件部分解释

Example

# 摘自李老师的gmx_MMPBSA.bsh程序

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
read 
    mol pqr pro.pqr
    mol pqr lig.pqr
    mol pqr com.pqr
end

# CALCULATE POTENTIAL FOR COMPLEX
elec name com
  mg-auto
  mol 3
  dime   193  161  161        # 格点数目, 所需内存: 954.199 MB
  cglen  245.472  215.334  224.289        # 粗略格点长度
  fglen  91.824  81.778  84.763        # 细密格点长度
  fgcent 41.621  32.573  46.6415  # 细密格点中心
  cgcent 41.621  32.573  46.6415  # 粗略格点中心

  temp  298.15      # 温度
  pdie  2           # 溶质介电常数
  sdie  78.54       # 溶剂介电常数, 真空1, 水78.54

  npbe              # PB方程求解方法, lpbe(线性), npbe(非线性), smbpe(大小修正)
  bcfl  mdh         # 粗略格点PB方程的边界条件, zero, sdh/mdh(single/multiple Debye-Huckel), focus, map
  srfm  smol        # 构建介质和离子边界的模型, mol(分子表面), smol(平滑分子表面), spl2/4(三次样条/7阶多项式)
  chgm  spl4        # 电荷映射到格点的方法, spl0/2/4, 三线性插值, 立方/四次B样条离散
  swin  0.3         # 立方样条的窗口值, 仅用于 srfm=spl2/4

  srad  1.4         # 溶剂探测半径
  sdens 10          # 表面密度, 每A^2的格点数, (srad=0)或(srfm=spl2/4)时不使用

  ion charge  1 conc 0.15 radius 0.95  # 阳离子的电荷, 浓度, 半径
  ion charge -1 conc 0.15 radius 1.81  # 阴离子

  calcforce  no
  calcenergy comps
end

APOLAR name comSAS
  mol 3

  temp  298.15 # 温度
  srfm  sacc   # 构建溶剂相关表面或体积的模型
  swin  0.3    # 立方样条窗口(A), 用于定义样条表面

  # SASA
  srad  1.4    # 探测半径(A)
  gamma 1      # 表面张力(kJ/mol-A^2)

  #gamma const 0.027     0        # 表面张力, 常数
  #gamma const 0.0226778 3.84928  # 表面张力, 常数

  press  0     # 压力(kJ/mol-A^3)
  bconc  0     # 溶剂本体密度(A^3)
  sdens 10
  dpos  0.2
  grid  0.1 0.1 0.1

  # SAV
  #srad  1.29      # SAV探测半径(A)
  #press 0.234304  # 压力(kJ/mol-A^3)

  # WCA
  #srad   1.25           # 探测半径(A)
  #sdens  200            # 表面的格点密度(1/A)
  #dpos   0.05           # 表面积导数的计算步长
  #bconc  0.033428       # 溶剂本体密度(A^3)
  #grid   0.45 0.45 0.45 # 算体积分时的格点间距(A)

  calcforce no
  calcenergy total
end

print elecEnergy com end
print apolEnergy comSAS end

在此仅说明旧的APBS输入文档READ input file section — APBS 3.4.1 documentation

输入文件一共分为四个部分：

  • READ
  • ELEC
  • APOLAR
  • PRINT

READ

通用格式：
READ
[keywords…]
END

读入的分子按顺序映射id(1,2,3…)
换句话说，每个输入的pqr文件按照顺序划分id，从1开始，递增1

mol

读取分子数据 (最常用)

1
2
3
READ
	mol {format} {path}
END
  • format
    • pqr
    • pdb 需要指定parm参数体，提供原子的电荷和半径
  • path 路径

parm

指定原子电荷和半径，与mol读取pdb连用

1
2
3
READ
	parm {format} {path} 
END
  • format
    • flat 以APBS flat 格式指定
    • xml 以APBS XML 格式
  • path

charge

此命令允许apbs读取映射到网格的分子的固定电荷密度函数。输入为电荷密度的map。这些值的单位为$e_{c}\dot{A}^{-3}$,$e_{c}$为电子电荷常数。该命令将会读取由ELEC输出的电荷密度图。

1
2
3
READ 
		charge {format} {path} 
END
  • format 指定电荷图格式，允许的格式如下
    • dx OpenDX scalar data format
    • gz gzipped (zlib) compressed OpenDX scalar data format. (即可以直接读取压缩文件)
  • path 路径

diel

This command allows APBS to read the dielectric function mapped to 3 meshes shifted by one-half grid spacing in the x, y, and z directions. The inputs are maps of dielectric variables between the solvent and biomolecular dielectric constants; these values are unitless. In general, this command will read dielectric maps written by by ELEC input file section write commands. The syntax of this command is:

1
2
3
READ
	diel {format} {path-x} {path-y} {path-z}
END
  • format
    • dx
    • gz
  • path-x x方向移动介电图位置
  • path-y
  • path-z

kappa

This command allows APBS to read the ion-accessibility function mapped to a mesh. The inputs are maps of ion accessibility values which range between 0 and the build Debye-Hückel screening parameter; these values have units of Å-2. In general, this command will read kappa-maps written by by ELEC input file section write commands. The syntax of this command is:

1
2
3
READ
	kappa {format} {path}
END
  • format
    • dx
    • gz
  • path

此命令与diel连用

pot

读取静电势，输入的通常是先前计算中保存的静电势

1
2
3
READ
	pot {format} {path} 
END
  • format
    • dx
    • gz
  • path

ELEC

ELEC代码块用于计算极性溶剂化能(静电)计算

1
2
3
4
ELEC [ name {id} ]
	{type}
	{keyword...}
END

{id}变量指定当前计算任务名称(自定义)。当计算的任务过多时，指定名称是必要的。

{type}命令定义了ELEC计算的类型：

  • 使用PMG进行有限差分网格计算 (常用)
    • mg-auto
    • mg-para
    • mg-manual
  • 基于Geometric flow solvation的有限差分计算
    • geoflow-auto
  • 使用TABI-PB的边界元素方法计算
    • tabi
  • 解析和半解析解的PB近似
    • pbam-auto
    • pbsam-auto
  • FET有限元计算
    • fe-manual
  • No-op modes for generating coefficient maps
    • mg-dummy

{keyword}进一步指定特定类型的计算命令。在此指定极性溶剂化自由能计算的主要组成。in文件中可能会有很多个ELEC代码块，用于使用不同参数计算不同分子。某些类型的边界条件(bcfl)需要先前计算的信息，因此ELEC语句的顺序十分重要。

mg-auto

自动配置有限差分泊松玻尔兹曼计算

该计算自动设置并执行一系列的单点PBE计算，以聚焦于系统某一区域(结合位点)的计算。该命令可以被认为是mg-manual的自动执行版本。推荐使用此type

该方法先以粗糙网格求解PB方程，解被用于设置小区域的dirichlet边界条件值。然后使用一个finer格点计算感兴趣区域的其他部分。

除非另作说明，否则必须包含所有的关键字

keyword

mol

指定对哪个分子进行计算
mol {id}

bcfl

指定求解PB方程的边界条件类型. Usage: bcfl {flag}

  • zero
    边界势能设为0, 此条件不经常使用，如果不合理使用会导致巨大的错误
  • sdh
    “Single Debye-Hückel” boundary condition.
    边界势能设置为Debye-Huckel模型(一个带有点电荷，偶极子和四极得球体)预描述的值
    模型得球体半径被设置为生物分子的半径，球体电荷，偶极矩以及四极矩被设为蛋白的总和。
    当边界离生物分子足够远时，此边界条件最有效
  • mdh
    “Multiple Debye-Hückel” boundary condition
    边界势能设置为Debye-Huckel模型(带有点电荷的多个非相互作用球)预描述的值
    非接触球体的半径被设置为原子半径，球体电荷设置为原子电荷。
    当边界离生物分子较近时，此边界条件比sdh有效，对于生物大分子来说速度可能较慢
  • focus
    “Focusing” boundary condition.
    边界势能由先前的PB计算得出。
    此选项只在mg-manual类型中使用
    所有的边界点均在先前计算的区域内。如果任意边界点处于区域之外，将使用Debye-Huckel边界条件计算其值。
  • map
    Specifying map allows a previously calculated potential map to be used in a new focusing calculation. A typical scenario is using the same coarse grid for multiple focusing calculations. A potential map can be written once from a coarse grid calculation, then used in subsequent runs to bypass the need to recalculate the coarse grid
calcenergy

控制输出非极性溶剂自由能计算
calcenergy <flag>

  • no
    不计算任何能量(废弃)
  • total
    计算并返回整个分子的非极性能量
  • comps
    计算并返回整个分子的非极性能量以及每个原子的能量贡献
calcforce

控制非极性溶剂自由能能得输出
calcforce <flag>

  • no
    不计算任何力(废废弃)
  • total
    计算并返回整个分子的静电势和非极性力
  • comps
    计算并返回整个分子的静电势和非极性力以及每个原子的力贡献

以下为输出控制

  • tot {n}
    原子n的力
  • qf {n}
    原子n的固定电荷力
  • db {n}
    原子n的边界介电力
  • ib {n}
    原子n的离子边界力

输出的值将会以三列的形式打印，其对应于x，y，z方向的力矢量

cgcent

控制PB方程计算的静电势能输出
cgcent {mol id | xcent ycent zcent}

  • mol id
    将指定id的分子至于网格中心
  • xcent ycent zcent
    将网格居中在指定位置(xcent ycent zcent) (A), 基于PDB坐标帧
cglen

指定自动多网格PB计算的粗糙网格长度(mg-auto, mg-para)，在每个方向上可能有所不同
cglen {xlen ylen zlen}
该命令为起始网格设置，因此必须确保网格足够大，以能够包裹整个分子，并确保所选边界条件合适。

  • xlen ylen zlen

    在x, y ,z方向上的网格长度。(浮点数，A)

chgm

指定生物大分子电荷计算方法(mg-manual,mg-auto,mg-para)
chgm {flag}

  • sp10
    传统的三点插值法(线性样条曲线)。电荷被映射到最近邻的网格格点上，产生的电势对于网格间距、长度和位置非常敏感
  • sp12
    立方B样条曲线离散化。电荷被映射到最近和下一个最近的网格格点上。对于网格间距、长度和位置敏感性较小
  • sp14
    五次B样条曲线离散化。与sp12相似，除了电荷、多重偶极子也被映射到了下下个近邻网格格点上(125个格点接受电荷密度)
dime

指定基于网格的离散化每个核处理的网格点数

dime {nx ny nz}
对于mg-manual计算，参数取决于nlev的选择：$n=c2^{l+1}+1$, n为dime参数，c为非零整数， l为nlev值。
最常用的值为65，97，129，161(在不同的方向上可以指定不同的值)。
这些值均与nlev=4兼容
如果指定了不合适的value，APBS将会自动向下调整dime值到一个合适位置。这意味着不合适的value不会报错，但是会带来不精确的结果。

  • nx ny nz
    x，y，z方向上的格点数(integer)
etol

指定偏微分方程求解器迭代的公差
etol {tol}

  • tol
    误差限度(float)

该关键字可选，主要针对mg-manual, mg-auto, mg-para

fgcent

基于分子中心或者绝对坐标指定细网格中心 （mg-para, mg-auto）

fgcent {mol {id} | xcent ycent zcent}

同cgcent定义

  • mol {id}
  • xcent ycent zcent
fglen

指定自动多网格PB计算的细网格长度(mg-auto, mg-para)，在每个方向上可能有所不同
fglen {xlen ylen zlen}
该命令为起始网格设置，因此必须确保网格足够大，以能够包裹整个分子，并确保所选边界条件合适。

  • xlen ylen zlen

    在x, y ,z方向上的网格长度。(浮点数，A)

ion

指定体系中存在的离子浓度
可重复此命令指定多种类型的离子 ，但是仅使用最大的离子半径来决定离子可及性函数。
整个体系的离子必须呈现电中性，即设定的阴阳离子电荷密度和浓度必须相等
ion charge {charge} conc {conc} radius {radius}

  • charge

    离子电荷(浮点数)

  • conc
    离子浓度(浮点数，M)

  • radius
    离子半径(浮点数，A)

lpbe

指定求解线性泊松玻尔兹曼方程

lrpbe

指定求解正则化泊松玻尔兹曼方程的线性格式。
正则化的PBE方程替代描述点电荷分布的格林方程。
求解结果对应反应场而不是整体势能。整体势能通过添加一定的库伦项获得。
该方程可以立即产生溶剂化能，而无需参考计算。

npbe

求解非线性泊松玻尔兹曼方程

pdie

指定溶质分子的介电常数
pdie {diel}

diel为无单位的生物分子介电常数(float).
取值为2-20。取值较低时仅考虑电子极化，取值较高时同时考虑因为分子内运动导致的额外极化。
介电常数必须大于1.

sdens

指定计算表面项(分子表面、溶剂可及表面积)时每A2范围内的格点数。
当srad=0.0或者srfm=spl2时，自动忽略

sdens {density}

默认值为10.0

sdie

指定溶剂的介电常数
sdie {diel}
diel用于描述溶剂的介电常数(浮点数)。
数值必须大于1，生物体系下的水的介电常数通常在78-80之间。

srad

指定溶剂分子半径
该参数用于定义与溶剂相关的表面和体积。
对于水分子的分子表面来说，设置为1.4 A
对于vdw表面来说，设置为0 A
srad {radius}

在命令参数srfm spl2存在下，该选项失效

srfm (elec)

指定用于构建介电常数和离子可及性系数的模型
srfm {flag}

  • mol
    介电系数根据分子表面并以而定义
  • smol
    同上，另外，通过9点谐波平均来平滑，以期在某种程度上降低对网格设置的敏感性 (default)
  • spl2
    介电常数和离子可及行系数由立方样条曲线表面定义
    介电窗口由swim控制
  • spl4
    介电常数和离子可及行系数由七阶多项式定义
swin

立方样条的窗口值
swin {win}
win为浮点数，样条窗口宽度 通常设置为0.3A

如果更改了样条窗口，通常需要修改力场参数

temp

指定计算温度
temp {t}
浮点数，单位 K

usemap

指定PB计算中应用的预计算的系数图，必须先在READ代码块中指定。
usemap {type} {id}

  • type
    • diel
      Dielectric function map (as read by READ input file section diel); this causes the pdie, sdie, srad, swin, and srfm (elec) parameters and the radii of the biomolecular atoms to be ignored when computing dielectric maps for the Poisson-Boltzmann equation. Note that the pdie and sdie values are still used for some boundary condition calculations as specified by bcfl.
    • kappa
      Mobile ion-accessibility function map (as read by READ input file section kappa); this causes the swin and srfm (elec) parameters and the radii of the biomolecular atoms to be ignored when computing mobile ion values for the Poisson-Boltzmann equation. The ion parameter is not ignored and will still be used.
    • charge
      Charge distribution map (as read by READ input file section charge); this causes the chgm parameter and the charges of the biomolecular atoms to be ignored when assembling the fixed charge distribution for the Poisson-Boltzmann equation.
    • pot
      Potential map (as read by READ input file section pot); this option requires setting bcfl to map.
  • id
write

控制PB运算期间标量数据的输出。
可以重复多次，以提供来自APBS各种类型数据的输出
write {type} {format} {stem}

  • type

    • charge
      输出在每A3空间内生物分子的电荷分布

    • pot
      输出整个problem域内的静电势能，单位kb T ec-1 
      kb
      玻尔兹曼常数：1.3806504 x 10-23 J K-1

      T
      设定温度

      ec-1
      电子电荷 1.60217646 × 10-19 C

    • atompot
      每个原子位置的静电势

    • smol
      输出由分子表面定义的溶剂可及性，与 srfm (elec) smol连用
      值从0(不可及)到1(可及)

    • sspl
      输出基于样条曲线的溶剂可及性，与 srfm (elec) spl2连用

    • vdw
      输出基于范德华的溶剂可及性，与 srfm (elec) smol和srad 0.0连用

    • ivdw
      输出基于范德华的溶剂可及性，与 srfm (elec) smol连用

    • lap
      Write out the Laplacian of the potential $∇^2ϕ$ in units of $k_BTe_c^{-1} Å^{-2}$

    • edens
      Write out the “energy density” −∇⋅ϵ∇ϕ−∇⋅ϵ∇ϕ in units of $k_BTe_c^{-1} Å^{-2}$

    • ndens
      输出离子数值密度

    • qdens
      输出离子电荷密度

    • dielx | diely | dielz
      Write out the dielectric map shifted by 1/2 grid spacing in the {x, y, z}-direction

  • format

    • dx
      Write out data in OpenDX scalar data format.
    • avs
      Write out data in AVS UCD format.
    • uhbd
      Write out data in UHBD scalar data format.
    • gz
      Write out OpenDX scalar data format in gzipped (zlib) compatible format. Appends .dx.gz to the filename.
    • flat
      Write out data as a plain text file

  • stem
    输出路径

writemat

控制PB方程的数学运算符输出(基于Harwell-Boeing矩阵格式)
writemat {type} {stem}

  • type
    • poisson
      输出泊松算符$-\nabla \cdot \epsilon \nabla$
  • stem
    输出路径

mg-para

其他选项同mg-auto

async

执行异步并行
async {rank}

ofrac

指定并行运算中单核网格重叠量
ofrac {frac}
where frac is a floating point value between 0.0 and 1.0 denoting the amount of overlap between processors. Empirical evidence suggests that an value of 0.1 is sufficient to generate stable energies.

pdime

指定并行focus运算时的处理器数组
pdime {npx npy npz}

mg-manual

gcent

基于分子中心或者绝对坐标指定网格中心 **mg-manual **
fgcent {mol id | xcent ycent zcent}

  • mol {id}
  • xcent ycent zcent
glen

指定网格域长度 **mg-manual **
glen {xlen ylen zlen}

grid

指定网格间隔 **mg-manual **
grid hx hy hz

nlev

指定多级层次结构的深度 **mg-manual **

nlev {lev}

geoflow-auto

tabi

pbam-auto

pbsam-auto

fe-manual

mg-dummy

APOLAR

此代码块用于计算非极化溶剂自由能。

1
2
3
APOLAR [name id]
	<keywords...>
END

非极性溶剂化自由能有下述公式计算：
$$
W^{np}(x)=\gamma A(x)+pV(X)+\overline{\rho}\sum_{i=1}^{N}\int_{\Omega}{u_{i}^{att}(x_{i},y)\theta)(x,y)}dy
$$
平均非极性溶剂化力：
$$
F_{i}^{(np)}(x)=-\gamma \frac{\partial A(x)}{\partial{x_{i}}}-p\int_{\Gamma (x)}{\frac{y-x_{i}}{||y-x_{i}||}}dy-\overline{\rho}\sum_{i=1}^{N}\int_{\Omega}\frac{\partial u_{i}^{(att)}(x_{i},y)}{\partial x_{i}}\theta (x,y)dy
$$
$\gamma$：溶剂排斥表面张力(硬球模型) (gamma)

$A$：与构象相关的溶质表面积 (srad srfm(apolar))

$p$：排斥溶剂压力(硬球模型) (press)

$V$: 与构象有关的溶质体积 (srad srfm(apolar))

$\rho$：溶剂密度 (bconc)

bconc

指定溶剂密度。

设置为0，即消去非极性溶剂化计算的积分贡献

bconc <density>

单位：$A^{-3}$

calcenergy

控制非极性溶剂化自由能计算的能量输出

calcenergy <flag>

  • no

    不计算能量

  • total

    计算并返回整体非极性能量

  • comps

    计算并返回整体非极性能量和每个原子的能量贡献

calcforce

控制非极性溶剂自由能能得输出
calcforce <flag>

  • no
    不计算任何力(废废弃)
  • total
    计算并返回整个分子的静电势和非极性力
  • comps
    计算并返回整个分子的静电势和非极性力以及每个原子的力贡献

以下为输出控制

  • tot {n}
    原子n的力
  • qf {n}
    原子n的固定电荷力
  • db {n}
    原子n的边界介电力
  • ib {n}
    原子n的离子边界力

输出的值将会以三列的形式打印，其对应于x，y，z方向的力矢量

dpos

用于表面积有限差分计算的位移

dpos {displacement}

gamma

指定表面张力系数

gamma <value>

grid

指定网格间隔
grid hx hy hz

mol

指定对哪个分子进行计算
mol {id}

press

溶剂压力 $ kJ mol^{-1} Å^{-3}$

press {value}

sdens

指定计算表面项(分子表面、溶剂可及表面积)时每A2范围内的格点数。
当srad=0.0或者srfm=spl2时，自动忽略

sdens {density}

默认值为10.0

srad

指定溶剂分子半径
该参数用于定义与溶剂相关的表面和体积。
对于水分子的分子表面来说，设置为1.4 A
对于vdw表面来说，设置为0 A
srad {radius}

在命令参数srfm spl2存在下，该选项失效

srfm(apolar)

指定用于构建介电常数和离子可及性系数的模型
srfm {flag}

  • sacc
    溶剂可及表面和体积

swin

立方样条的窗口值
swin {win}
win为浮点数，样条窗口宽度 通常设置为0.3A

如果更改了样条窗口，通常需要修改力场参数

temp

指定计算温度
temp {t}
浮点数，单位 K

PRINT

该section允许将计算结果进行线性组合并标准输出

1
PRINT {what} [id op id op ...] END

{what}是强制指定的：

  • elecEnergy 输出先前ELEC计算的静电势能
  • elecForce 输出先前ELEC计算的静电力
  • apolEnergy 输出先前APOLAR计算的非极性能量
  • apolForce 输出先前APOLAR计算的非极性力

紧接着是一系列的id op id op id op...，id即为之前指定name

  • id 字符串/数字，指定ELEC或者APOLAR代码块的id
  • op 算术操作符， + ，-
1
2
3
4
5
6
7
# Example
# Energy change due to binding
print energy complex - ligand - protein end
# Energy change due to solvation
print energy solvated - reference end
# Solvation energy change due to binding
print energy complex_solv - complex_ref - ligand_solv + ligand_ref - protein_solv + protein_ref end