【lammps】固液界面の密度分布計算について

## 固液界面の密度分布
固体上に置かれた液体は界面の近くで密度分布に変化があるらしい。

- 固体-液体界面での水分子分布の直接可視化に成功[京都大学]https://www.kyoto-u.ac.jp/sites/default/files/embed/jaresearchresearch_results2017documents171213_101.pdf
- 液中AFMによるマイカ基板表面の水和構造観察[材料科学技術振興財団]https://www.mst.or.jp/casestudy/tabid/1318/pdid/574/Default.aspx
- サファイア(001)と α-クォーツ(100)の表面構造および水和構造を原子スケールで可視化！[金沢大学]https://www.kanazawa-u.ac.jp/rd/129036

## LAMMPSでシミュレーションしてみる
lammpsを使って、固体と液体を作っていく。
固体は２次元に並べた原子の壁、液体はメタノールとする。
```
variable nstep index 100000

#--------------------------------基本設定----------------------------------
units real
atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style opls
improper_style cvff
pair_style lj/cut/coul/long 12 12
pair_modify mix arithmetic
special_bonds lj/coul 0.0 0.0 0.5
kspace_style pppm 0.0001


#read_data "methanol.settings"	

variable N_A  index 6.022140857e+23
variable mass_C index 12.001
variable mass_O index 15.9994
variable mass_H index 1.008


~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
分子の配置方法は省略

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

#------------------------メタノールの力場パラメータを設定------------------------------


#Pair	Coeffs				
					
pair_coeff 1 1 0.066 3.5	
pair_coeff 2 2 0.17 3.12	
pair_coeff 3 3 0.03	2.5	
pair_coeff 4 4 0.03	2.5	
pair_coeff 5 5 0.03 2.5	
pair_coeff 6 6 0 0	
					
#Bond	Coeffs				
					
bond_coeff 1 320 1.41		
bond_coeff 2 340 1.09		
bond_coeff 3 340 1.09		
bond_coeff 4 340 1.09		
bond_coeff 5 553 0.945		
					
#Angle	Coeffs				
					
angle_coeff	1	35	109.5		
angle_coeff	2	35	109.5		
angle_coeff	3	35	109.5		
angle_coeff	4	55	108.5		
angle_coeff	5	33	107.8		
angle_coeff	6	33	107.8		
angle_coeff	7	33	107.8		
					
#Dihedral	Coeffs				
					
dihedral_coeff	1	0	0	0.352	0
dihedral_coeff	2	0	0	0.352	0
dihedral_coeff	3	0	0	0.352	0
					
#Improper	Coeffs				
					
improper_coeff	1	0	-1	2	
improper_coeff	2	0	-1	2

#Masses

mass    1     12.011  
mass     2     15.999  
mass     3      1.008  
mass     4      1.008  
mass     5      1.008  
mass     6      1.008

#------------------------壁の力場パラメータを設定------------------------------
#Pair_style      lj/cut 
pair_coeff      7 7 1 1

#Masses
mass            7 1.0

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
分子の配置方法は省略

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

#---------------------------可視化用データの設定--------------------------------
compute Dens_liquid2D liquid chunk/atom bin/2d z lower 0.5 y lower 0.5
compute Dens_liquid liquid chunk/atom bin/1d z lower 0.005 units reduced
fix 2 liquid  ave/chunk  1 ${nstep} ${nstep} Dens_liquid density/mass ave running file liquid_density.dat
fix 3 liquid ave/chunk   1 ${nstep} ${nstep} Dens_liquid2D density/mass ave running file liquid_density2D.dat


#--------------------------------計算実行--------------------------------------

fix 1 system nve
fix vliquid liquid temp/rescale 10 300 300 0 1
compute comliquid liquid temp
dump methDump all atom 1000 dump.lammpstrj

timestep 1

thermo_style    custom step time temp press vol density c_comliquid #v_GamNsurf
thermo 1000

run ${nstep}
write_data data.last

```

 ## 計算結果まとめ
 ### 1次元データの解析
 以下の図はz方向でのメタノール密度分布を表している。図より、界面付近に最大の密度ピークが存在していることがわかる。この密度ピークは強度が小さくなるものの、周期的に発生していることが確認できる。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2371660/99708178-5502-c373-34c2-01063e09d492.png)

### 2次元データの解析
2次元データからpythonを用いて可視化した結果、1次元データでも見られた周期的な高密度領域が確認できた。また、界面に対して水平方向には密度に変化がないことがこの図からわかった。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2371660/b5470bf3-2e61-f8be-aedc-2e5b0d8571f7.png)

この周期的に高密度領域が現れる現象は実験でも観察が行われている。


## 参考文献
- https://pubs.aip.org/aip/jcp/article-abstract/138/6/064905/192609/Parameter-passing-between-molecular-dynamics-and?redirectedFrom=fulltext