martini粗粒化聚合物模型在lammps中使用fix bond/create指令聚合成键总是失败怎么办？

anananante

求助各位大佬，第一次做粗粒化，选用了martini力场，需要异氰酸酯和醇聚合。
我的键结参数都是martini官方教程中给出的bartender拟合生成的，pair coeffs遵循了martini3的粗粒化规则，目前就是单体盒子事先跑npt和nvt没有问题，一旦开始用fix bond/create指令成键，就会在运行一段时间后报错bond atom missing。我的体系按照数量最多能成2000个键，我调整过各种参数最多跑到成三四百个键左右就会报错。
请问是我的模型构建的有问题还是参数设置不合理？还是说martini有什么特殊的设置要求？下面贴上我的映射和in文件
目前是仅仅在令IPDI和PPG反应，OPPEA准备二阶反应，暂时还没用到，盒子中有这3个单体
# ----------------- Init Section -----------------
units           real
dimension 3
boundary        p p p
neighbor 3 bin
neigh_modify every 1 delay 0 check yes

#-----style-----------
atom_style      full
pair_style      lj/cut  20.0
pair_modify             tail no
pair_modify     mix arithmetic
bond_style      harmonic
angle_style     harmonic
special_bonds lj/coul 0.0 0.0 1.0 angle yes dihedral yes

# ----------------- Atom Definition Section -----------------
read_data stage1.data
# ============= coeffes Section ===================
pair_coeff     1        1                     0.88191    4.7
pair_coeff     2        2                     0.88191    4.7
pair_coeff     3        3                     0.81021    4.7
pair_coeff     4        4                     0.56165    4.1
pair_coeff     5        5                     0.67876    4.1
pair_coeff     6        6                0.6214     4.1
pair_coeff     7        7                     0.88191    4.7
pair_coeff     8        8                     0.36089    3.4
pair_coeff     9        9                     0.88191    4.7
pair_coeff     10        10                     0.67876    4.1

# -- harmonic bondType parameter  --
bond_coeff      1                 7.3439    3.45
bond_coeff      2                 7.3439    3.45
bond_coeff      3                 7.3439    3.45
bond_coeff      4                 7.3439    3.45
bond_coeff      5                 7.3439    3.45
bond_coeff      6                 7.3439    3.45
bond_coeff      7                 7.3439    3.45
bond_coeff      8                 7.3439    3.45
bond_coeff      9                 7.3439    3.45

# -- cosine/squared angleType parameter  --
angle_coeff     1                5.6804  180
angle_coeff     2                5.6804  180
angle_coeff     3                5.6804  180
angle_coeff     4                5.6804  180
angle_coeff     5                5.6804  180
angle_coeff     6                5.6804  180
angle_coeff     7                5.6804  180
angle_coeff     8                5.6804  180

# -- harmonic bondType parameter  --
bond_coeff      10                 62.9010      3.24
bond_coeff      11                 111.7513     2.55
bond_coeff      12                 4.6309       3.08
bond_coeff      13                 85.6537      2.77
bond_coeff      14                 125.9143     2.81

# -- cosine/squared angleType parameter  --
angle_coeff     9                3.2886   125.4
angle_coeff     10                7.9372   115.28
angle_coeff     11                4.1634   119.19
angle_coeff     12                8.3339   109.88


# -- harmonic bondType parameter  --
bond_coeff      15                 6.6422    3.41
bond_coeff      16                 4.5462    3.07
bond_coeff      17                 282.3161  2.08
bond_coeff      18                 294.4805  2.09
bond_coeff      19                 33.2672   2.36
bond_coeff      20                 296.0863  2.07
bond_coeff      21                 269.0198  2.1
bond_coeff      22                 273.3664  2.1
bond_coeff      23                 292.1409  2.07
bond_coeff      24                 7.0162    3.2
bond_coeff      25                 7.0162    3.04
# -- cosine/squared angleType parameter  --
angle_coeff     13                5.5137   125.26
angle_coeff     14                12.8128  103.29
angle_coeff     15                23.0826   92.45
angle_coeff     16                27.9845  145.31
angle_coeff     17                34.0288  106.72
angle_coeff     18                35.1402  163.08
angle_coeff     19                36.0842  161.65
angle_coeff     20                233.5891  104.79
angle_coeff     21                13.7640    180.0
angle_coeff     22                13.7640    180.0

#-----NVT(298K)--------------
#velocity all create 298 161013 rot yes dist gaussian
min_style sd
minimize 1.0e-4 1.0e-6  10000 1000000
timestep 1

fix stage1h all  bond/create 200  1 5  4.2  24  prob 0.05 16101 iparam 1 9  jparam  1 10 #atype 21
dump    mydump all  custom 10000 stage1.lammpstrj id type x y z mol
fix 1 all nvt temp 298 298 100  
restart 10000 tmp.restart

thermo_style custom step elapsed temp press pe ke etotal vol density f_stage1h[1]  f_stage1h[2]
thermo 1000

run  10000000
write_data stage1out.data