lammps建模——EMC

用ms建模经常遇到些脑溢血的问题：AC模块完全跑不动，还有就是msi2lmp这个程序总是会丢势函数的参数。最后只好换个方法用了EMC，20MB的大小在建模这块简直吊打MS。

1.EMC下载安装

EMC文件下载地址：python调用EMC程序生成lammps的data文件。-EMC文档类资源-CSDN下载（需要安装python）

下载完解压后，在终端输入以下命令安装调用EMC的python库：

pip install emc-pypi

在python中执行以下代码，没有报错则说明改库安装成功。

import pyemc

由于EMC主体是用perl写的，所以电脑还需要支持perl环境，在终端输入以下命令查看是否支持perl：

perl -v

输出如下则表明已安装perl：

如果提示找不到perl，自行上官网下载安装就行了。

2.esh文件编写

下载的文件中自带了一些esh的例子，照着修改可以了。我用于生成PET/PE复合材料的esh如下：

#!/usr/bin/env emc_setup.pl

ITEM    OPTIONS

replace    true #覆盖存在的脚本
mass       true #在data中添加质量
ntotal     7200 #体系原子总数
density    0.3  #体系密度
field      pcff #设置力场（EMC支持COMPASS, PCFF, CHARMM, OPLS, TraPPE or coarse-grained force fields）
build_dir  .    #我也不知道有啥用

ITEM    END

ITEM    SHORTHAND #提供了一种快速构建化学文件的方法

PET O=C(C1=CC=C(C(OCCOC(C2=CC=C(C(OCCOC(C3=CC=C(C(OCCOC(C4=CC=C(C(OCCOC(C5=CC=C(C(OCCOC(C6=CC=C(C(OCCOC(C7=CC=C(C(OCCOC(C8=CC=C(C(OCCOC(C9=CC=C(C(OCCOC(C%10=CC=C(C(C)=O)C=C%10)=O)=O)C=C9)=O)=O)C=C8)=O)=O)C=C7)=O)=O)C=C6)=O)=O)C=C5)=O)=O)C=C4)=O)=O)C=C3)=O)=O)C=C2)=O)=O)C=C1)OCCOC,10 #名称，SMILE式，链数
PE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC,20 #名称，SMILE式，链数

ITEM    END

ITEM    COMMENTS #设置体系组分，含量

ITEM    GROUPS #组分
PET   O=C(C1=CC=C(C(OCCOC(C2=CC=C(C(OCCOC(C3=CC=C(C(OCCOC(C4=CC=C(C(OCCOC(C5=CC=C(C(OCCOC(C6=CC=C(C(OCCOC(C7=CC=C(C(OCCOC(C8=CC=C(C(OCCOC(C9=CC=C(C(OCCOC(C%10=CC=C(C(C)=O)C=C%10)=O)=O)C=C9)=O)=O)C=C8)=O)=O)C=C7)=O)=O)C=C6)=O)=O)C=C5)=O)=O)C=C4)=O)=O)C=C3)=O)=O)C=C2)=O)=O)C=C1)OCCOC #名称，smile式
PE CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC #名称，smile式
ITEM    END

ITEM    CLUSTERS #含量

PET   PET,5 #名称，名称，链数
PE PE,20 #名称，名称，链数

ITEM    END #含量

将上面这段代码保存为1.esh。

smile式在许多化学资料库上可以查找，但是还是会有很多缺少的。建议还是自己用chemdraw绘制并导出smile式。

注释自己删掉

3.执行EMC生成data文件

在python中运行以下程序

import pyemc

pyemc.setup('xxx.esh')
pyemc.build('build.emc')

在esh文件所在目录下会生成以下文件：

4.文件详情

params文件中主要包含模型的势函数，质量等参数:

......
# Masses

mass		1    12.01115  # c
mass		2    12.01115  # c_0
mass		3    12.01115  # c_1
mass		4    12.01115  # cp
mass		5     1.00797  # hc
mass		6    15.99940  # o_1
mass		7    15.99940  # o_2
mass		8    15.99940  # oc

# Potentials

pair_style	lj/class2/coul/long ${cutoff} ${charge_cutoff}  # 9.5
bond_style	class2
angle_style	class2
dihedral_style	class2
improper_style	class2
pair_modify	mix sixthpower tail yes
special_bonds	lj/coul 0 0 1

# Pair Coeffs

pair_coeff	 1  1    0.05400    4.01000  # c,c
pair_coeff	 2  2    0.12000    3.30800  # c_0,c_0
pair_coeff	 3  3    0.12000    3.81000  # c_1,c_1
pair_coeff	 4  4    0.06400    4.01000  # cp,cp
pair_coeff	 5  5    0.02000    2.99500  # hc,hc
pair_coeff	 6  6    0.26700    3.30000  # o_1,o_1
pair_coeff	 7  7    0.24000    3.42000  # o_2,o_2
pair_coeff	 8  8    0.24000    3.53500  # oc,oc
......

data文件详情如下：

LAMMPS output created by EMC v9.4.4, build Aug 20 2021 07:52:52

        7140  atoms
        7185  bonds
       13760  angles
       19981  dihedrals
        8520  impropers

           8  atom types
          12  bond types
          20  angle types
          26  dihedral types
          15  improper types

             0    60.21913826  xlo xhi
             0    60.21913826  ylo yhi
             0    60.21913826  zlo zhi

Masses

       1   12.01115  # c
       2   12.01115  # c_0
       3   12.01115  # c_1
       4   12.01115  # cp
       5    1.00797  # hc
       6   15.99940  # o_1
       7   15.99940  # o_2
       8   15.99940  # oc

Atoms

       1       1   6 -0.5310  42.9929423661  30.7556164498   8.4423851246  # o_1
       2       1   3  0.7200  41.8686922004  31.1522224122   8.7061678569  # c_1
       3       1   4 -0.0180  40.5796544532  30.6967084270   8.0664396047  # cp
       4       1   4 -0.1268  40.6052652233  29.6647885167   7.1066772045  # cp
       5       1   5  0.1268  41.5743523233  29.1512023110   6.8297214122  # hc
......

将data导入ovito中：

2022.10.21更新:错误解决 

最近有几位同学来问的我程序执行报这样一个错误：

这个大概是github上的原作者更新pyemc库时写了啥bug，导致程序无法运行。

解决的办法就是安装老版本的pyemc库，然后用新版本包里的.py文件替换上面路径中提示的.py文件。