Dear Steve,
When I use the PPPM to simulate cubic CaF2 structure with Sigma5 interface, the atoms was lost without warning even though I use the thermo_modify lost warn.
Just before I use the nvt to increase the temperature from 1K to 300K, there are 15180 atoms (50600 Ca, 101200 F) and the charge of the system is neutral. After 200ps nvt relaxation, there are only 149746 atoms left and the system is no longer charge free. The volume of the structure and box lengths remain unchanged during the simulation. So I want to know what is wrong.
Log file:
…………………
group Ca type 1
50600 atoms in group Ca
group F type 2
101200 atoms in group F
set group Ca charge 2
50600 settings made for charge
set group F charge -1
101200 settings made for charge
mass 1 40.078
mass 2 18.9984
lattice none
lattice fcc $a orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
lattice fcc 5.43356 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
Lattice spacing in x,y,z = 5.43356 5.43356 5.43356
#Potential and atomic settings
kspace_style pppm 0.0001
pair_style buck/coul/long 10
pair_coeff 1 1 0.0 1.0 0.0
pair_coeff 1 2 674.3 0.336 0.0
pair_coeff 2 2 1808.0 0.293 109.1
neighbor 0.3 bin
neigh_modify delay 0 every 1 check yes
#Time step
timestep 0.001
Previous compute and variable
compute mytemp all temp/com
compute EKe all ke/atom
compute EPe all pe/atom
variable Them_MP_t equal 2000000*9/10
variable equl_t equal 2000000/10
variable Tprofile_t equal 2000000/20
variable resta_intevel equal 2000000/10
#constant parameters
variable kB equal 8.617343e-5
variable A2m equal 1.0e-10 # A --> m
variable ps2s equal 1.0e-12 # ps --> s
variable eV2J equal 1.60218e-19 # eV --> J
variable temp_sys equal 300
variable Nevery index 100
variable bin_MP equal ${slicetotal}
variable bin_MP equal 920
variable cross_sec equal lylz(v_A2m)^2
variable time equal (step-2*${equl_t}+10e-6)dtv_ps2s
variable time equal (step-2*200000+10e-6)dtv_ps2s
variable flux equal (f_4v_eV2J)/(2.0v_cross_sec*v_time)
compute strs all stress/atom virial
compute strs1 all stress/atom
compute Jflux1 J_region1 heat/flux EKe EPe strs # units: eV*(A/ps)
compute Jflux2 J_region2 heat/flux EKe EPe strs # units: eV*(A/ps)
variable subvol equal (bound(J_region1,xmax,region1)-bound(J_region1,xmin,region1))v_A2mv_cross_sec
print “${subvol}”
7.741943775e-25
variable J1x equal (c_Jflux1[1]v_eV2Jv_A2m/v_ps2s)/(v_subvol)
variable J1y equal (c_Jflux1[2]v_eV2Jv_A2m/v_ps2s)/(v_subvol)
variable J1z equal (c_Jflux1[3]v_eV2Jv_A2m/v_ps2s)/(v_subvol)
variable J2x equal (c_Jflux2[1]v_eV2Jv_A2m/v_ps2s)/(v_subvol)
variable J2y equal (c_Jflux2[2]v_eV2Jv_A2m/v_ps2s)/(v_subvol)
variable J2z equal (c_Jflux2[3]v_eV2Jv_A2m/v_ps2s)/(v_subvol)
variable J1 equal sqrt(v_J1xv_J1x+v_J1yv_J1y+v_J1z*v_J1z)
#temperature profile of selected sites
variable site_cold_l equal round({cold_l}/{slice_d})
variable site_cold_l equal round(96.85/${slice_d})
variable site_cold_l equal round(96.85/0.1118033989)
variable site_cold_r equal round({cold_r}/{slice_d})
variable site_cold_r equal round(6/${slice_d})
variable site_cold_r equal round(6/0.1118033989)
variable site_hot_l equal round({hot_l}/{slice_d})
variable site_hot_l equal round(45.425/${slice_d})
variable site_hot_l equal round(45.425/0.1118033989)
variable site_hot_r equal round({hot_r}/{slice_d})
variable site_hot_r equal round(57.425/${slice_d})
variable site_hot_r equal round(57.425/0.1118033989)
variable tCL equal f_30[${site_cold_l}][3]
variable tCL equal f_30[866][3]
variable pCL equal f_30[${site_cold_l}][4]
variable pCL equal f_30[866][4]
variable sCL equal f_30[${site_cold_l}][5]
variable sCL equal f_30[866][5]
variable tCR equal f_30[${site_cold_r}][3]
variable tCR equal f_30[54][3]
variable pCR equal f_30[${site_cold_r}][4]
variable pCR equal f_30[54][4]
variable sCR equal f_30[${site_cold_r}][5]
variable sCR equal f_30[54][5]
variable tHL equal f_30[${site_hot_l}][3]
variable tHL equal f_30[406][3]
variable pHL equal f_30[${site_hot_l}][4]
variable pHL equal f_30[406][4]
variable sHL equal f_30[${site_hot_l}][5]
variable sHL equal f_30[406][5]
variable tHR equal f_30[${site_hot_r}][3]
variable tHR equal f_30[514][3]
variable pHR equal f_30[${site_hot_r}][4]
variable pHR equal f_30[514][4]
variable sHR equal f_30[${site_hot_r}][5]
variable sHR equal f_30[514][5]
#atom’s temperature subtract the translation of the system
variable temp_atom atom (c_EKe-1/2mass(all)/count(all)(vcm(all,x)*vcm(all,x)+vcm(all,y)*vcm(all,y)+vcm(all,z)vcm(all,z)))/(1.5v_kB)
variable v equal vol
variable v0 equal $v
variable v0 equal 2029047.37
variable inix equal lx
variable lx0 equal ${inix}
variable lx0 equal 558.841646
variable iniy equal ly
variable ly0 equal ${iniy}
variable ly0 equal 60.7472008
variable iniz equal lz
variable lz0 equal ${iniz}
variable lz0 equal 59.76916
variable factor equal ${v0}/count(all)*10000
variable factor equal 2029047.37/count(all)*10000
variable Vatom equal ${v0}/count(all)
variable Vatom equal 2029047.37/count(all)
compute disp all displace/atom
compute ack all ackland/atom
compute RDF all rdf 100
compute diffu all msd
compute Energy all reduce ave c_EKe c_EPe
compute coord all coord/atom 2.3543
compute centro all centro/atom 4
compute cna all cna/atom 2.3543
compute V all reduce sum c_strs1[1] c_strs1[2] c_strs1[3] c_strs1[4] c_strs1[5] c_strs1[6] c_strs[1]
variable Vmis equal “sqrt(((c_V[1]-c_V[2])^2+(c_V[2]-c_V[3])^2+(c_V[3]-c_V[1])^2+6.0*((c_V[4])^2+(c_V[5])^2+(c_V[6])^2))/2.0)”
variable strss3 atom c_strs2[1]/${factor}
variable strss3 atom c_strs2[1]/133665.8347
variable s11 atom c_strs[1]/${factor}
variable s11 atom c_strs[1]/133665.8347
variable s22 atom c_strs[2]/${factor}
variable s22 atom c_strs[2]/133665.8347
variable s33 atom c_strs[3]/${factor}
variable s33 atom c_strs[3]/133665.8347
variable s12 atom c_strs[4]/${factor}
variable s12 atom c_strs[4]/133665.8347
variable s13 atom c_strs[5]/${factor}
variable s13 atom c_strs[5]/133665.8347
variable s23 atom c_strs[6]/${factor}
variable s23 atom c_strs[6]/133665.8347
variable mises atom “sqrt(((v_s11-v_s22)^2+(v_s22-v_s33)^2+(v_s33-v_s11)^2+6.0*((v_s12)^2+(v_s13)^2+(v_s23)^2))/2.0)”
Relaxation
thermo_style custom step temp etotal lx ly lz press vol
thermo_modify lost warn line one norm yes flush yes
thermo 5000
fix 100 all box/relax aniso 0.0 vmax 0.001
min_style cg
minimize 1.0e-5 1.0e-6 10000 5000
PPPM initialization …
G vector (1/distance) = 0.288635
grid = 288 50 50
stencil order = 5
RMS precision = 8.33632e-05
using single precision FFTs
brick FFT buffer size/proc = 69575 57600 27225
Memory usage per processor = 25.8541 Mbytes
Step Temp TotEng Lx Ly Lz Press Volume
0 0 -9.026047 558.84165 60.747201 59.76916 -5169.4169 2029047.4
15 0 -9.0280789 550.45902 60.939218 60.038044 -565.14698 2013948.7
Loop time of 7.04154 on 16 procs for 15 steps with 151800 atoms
Minimization stats:
Stopping criterion = energy tolerance
Energy initial, next-to-last, final =
-9.02604698029 -9.02799331782 -9.02807886119
Force two-norm initial, final = 28060.1 10766.8
Force max component initial, final = 27313.9 9303.21
Final line search alpha, max atom move = 9.63316e-08 0.000896193
Iterations, force evaluations = 15 15
Pair time (%) = 3.68537 (52.3376)
Kspce time (%) = 0.853661 (12.1232)
Neigh time (%) = 1.15532 (16.4072)
Comm time (%) = 0.0319651 (0.45395)
Outpt time (%) = 0 (0)
Other time (%) = 1.31522 (18.6781)
FFT time (% of Kspce) = 0.430923 (50.4794)
FFT Gflps 3d (1d only) = 6.76036 24.7671
Nlocal: 9487.5 ave 9570 max 9350 min
Histogram: 4 0 0 0 0 4 0 0 0 8
Nghost: 17943.1 ave 18092 max 17842 min
Histogram: 3 4 1 0 1 1 2 1 2 1
Neighs: 1.67277e+06 ave 1.69349e+06 max 1.64662e+06 min
Histogram: 4 0 0 0 4 0 0 2 4 2
FullNghs: 0 ave 0 max 0 min
Histogram: 16 0 0 0 0 0 0 0 0 0
Total # of neighbors = 26764320
Ave neighs/atom = 176.313
Neighbor list builds = 15
Dangerous builds = 15
unfix 100
reset_timestep 0
Equilibration
dump 1 all cfg 40000 tilt_CaF2_MP_snap.*.cfg id type xs ys zs vx vy vz v_temp_atom c_disp[4] c_ack c_centro c_cna c_EPe c_EKe v_mises
dump_modify 1 element Ca F
velocity all create 1.0 5812775 mom yes rot yes dist gaussian units box
fix 1 all nvt temp 1.0 ${temp_sys} 0.1 drag 2.0
fix 1 all nvt temp 1.0 300 0.1 drag 2.0
fix_modify 1 temp mytemp
fix 20 all ave/time 1 1000 1000 v_J2x v_J1x v_J1y v_J1z file thermo_CaF2_MP_equ1_${L_x}.txt
fix 20 all ave/time 1 1000 1000 v_J2x v_J1x v_J1y v_J1z file thermo_CaF2_MP_equ1_102.85.txt
#========thickness of the slice for calculat temperature profile=====================================
variable delta_s equal lx/a/{slicetotal}
variable delta_s equal lx/5.43356/${slicetotal}
variable delta_s equal lx/5.43356/920
fix 21 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower {delta_s} v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_{L_x}_${slice_d}.txt units lattice ave running
fix 21 all ave/spatial 1 100000 {Tprofile_t} x lower {delta_s} v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_{L_x}_{slice_d}.txt units lattice ave running
fix 21 all ave/spatial 1 100000 100000 x lower {delta_s} v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_{L_x}_${slice_d}.txt units lattice ave running
fix 21 all ave/spatial 1 100000 100000 x lower 0.1101165761 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_{L_x}_{slice_d}.txt units lattice ave running
fix 21 all ave/spatial 1 100000 100000 x lower 0.1101165761 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_102.85_${slice_d}.txt units lattice ave running
fix 21 all ave/spatial 1 100000 100000 x lower 0.1101165761 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_102.85_0.1118033989.txt units lattice ave running
#fix 22 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower 0.50 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_${L_x}_0.5.txt units lattice ave running
#fix 23 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower 0.75 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_${L_x}_0.75.txt units lattice ave running
#fix 24 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower 1.00 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_${L_x}_1.txt units lattice ave running
#fix 25 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower 2.00 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_${L_x}_2.txt units lattice ave running
#fix 26 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower 4.00 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_${L_x}_4.txt units lattice ave running
#fix 27 all ave/spatial 1 {Tprofile_t} {Tprofile_t} x lower 10.0 v_temp_atom c_EPe v_s11 file temp_profile_CaF2_MP_equ1_${L_x}_10.txt units lattice ave running
fix 30 all ave/spatial 1 1000 1000 x lower ${slice_d} v_temp_atom c_EPe v_s11 units lattice
fix 30 all ave/spatial 1 1000 1000 x lower 0.1118033989 v_temp_atom c_EPe v_s11 units lattice
fix 31 all ave/time 1 1 1000 v_tHL v_tHR v_tCL v_tCR v_pHL v_pHR v_pCL v_pCR v_sHL v_sHR v_sCL v_sCR file TPS_CaF2_MP_equ1_${L_x}.txt
fix 31 all ave/time 1 1 1000 v_tHL v_tHR v_tCL v_tCR v_pHL v_pHR v_pCL v_pCR v_sHL v_sHR v_sCL v_sCR file TPS_CaF2_MP_equ1_102.85.txt
restart {equl_t} CaF2_{L_x}_${temp_sys}.equil
restart 200000 CaF2_{L_x}_{temp_sys}.equil
restart 200000 CaF2_102.85_${temp_sys}.equil
restart 200000 CaF2_102.85_300.equil
run ${equl_t}
run 200000
PPPM initialization …
G vector (1/distance) = 0.287089
grid = 270 50 50
stencil order = 5
RMS precision = 9.14709e-05
using single precision FFTs
brick FFT buffer size/proc = 66550 54000 36300
Memory usage per processor = 68.3622 Mbytes
Step Temp TotEng Lx Ly Lz Press Volume
0 1 -9.0280566 550.45902 60.939218 60.038044 -615.82277 2013948.7
5000 9.0167807 -9.1007425 550.45902 60.939218 60.038044 12837.098 2013948.7
10000 15.496065 -9.0997859 550.45902 60.939218 60.038044 12926.743 2013948.7
15000 23.046913 -9.0982549 550.45902 60.939218 60.038044 13068.295 2013948.7
20000 30.547917 -9.0964225 550.45902 60.939218 60.038044 13259.67 2013948.7
25000 38.109393 -9.0944409 550.45902 60.939218 60.038044 13604.565 2013948.7
30000 45.485574 -9.0928581 550.45902 60.939218 60.038044 13815.193 2013948.7
35000 53.297713 -9.0906585 550.45902 60.939218 60.038044 14206.157 2013948.7
40000 60.161753 -9.0889569 550.45902 60.939218 60.038044 14404.328 2013948.7
45000 68.050519 -9.087208 550.45902 60.939218 60.038044 14590.241 2013948.7
50000 75.397999 -9.085399 550.45902 60.939218 60.038044 14924.506 2013948.7
55000 83.227336 -9.0832476 550.45902 60.939218 60.038044 15373.282 2013948.7
60000 90.037448 -9.0818082 550.45902 60.939218 60.038044 15350.132 2013948.7
65000 98.36613 -9.0797445 550.45902 60.939218 60.038044 15745.024 2013948.7
70000 104.91506 -9.0777649 550.45902 60.939218 60.038044 16170.024 2013948.7
75000 113.15339 -9.0760535 550.45902 60.939218 60.038044 16250.236 2013948.7
80000 120.19822 -9.0743654 550.45902 60.939218 60.038044 16531.282 2013948.7
85000 128.33161 -9.0724184 550.45902 60.939218 60.038044 16874.1 2013948.7
90000 134.76958 -9.0705461 550.45902 60.939218 60.038044 17117.597 2013948.7
95000 143.27988 -9.0689671 550.45902 60.939218 60.038044 17385.591 2013948.7
100000 149.86753 -9.0670376 550.45902 60.939218 60.038044 17636.338 2013948.7
105000 158.03651 -9.0651293 550.45902 60.939218 60.038044 17799.49 2013948.7
110000 165.25397 -9.063357 550.45902 60.939218 60.038044 18293.954 2013948.7
115000 173.32775 -9.0619549 550.45902 60.939218 60.038044 18372.202 2013948.7
120000 179.97175 -9.0594654 550.45902 60.939218 60.038044 18594.094 2013948.7
125000 187.3611 -9.0580838 550.45902 60.939218 60.038044 18948.314 2013948.7
130000 195.18919 -9.0564536 550.45902 60.939218 60.038044 19250.133 2013948.7
135000 203.11414 -9.0544041 550.45902 60.939218 60.038044 19277.163 2013948.7
140000 209.45954 -9.0521898 550.45902 60.939218 60.038044 19758.995 2013948.7
145000 217.56119 -9.051528 550.45902 60.939218 60.038044 19852.675 2013948.7
150000 224.6669 -9.0487551 550.45902 60.939218 60.038044 20190.172 2013948.7
155000 231.85752 -9.0469117 550.45902 60.939218 60.038044 20528.694 2013948.7
160000 240.35533 -9.0460201 550.45902 60.939218 60.038044 20551.791 2013948.7
165000 247.70477 -9.043755 550.45902 60.939218 60.038044 20852.143 2013948.7
170000 254.50526 -9.0410563 550.45902 60.939218 60.038044 21468.28 2013948.7
175000 262.1687 -9.0408383 550.45902 60.939218 60.038044 21294.243 2013948.7
180000 270.63241 -9.0382222 550.45902 60.939218 60.038044 21609.603 2013948.7
185000 276.26124 -9.0356376 550.45902 60.939218 60.038044 22193.569 2013948.7
190000 284.79471 -9.03511 550.45902 60.939218 60.038044 22127.274 2013948.7
195000 291.73344 -9.0331552 550.45902 60.939218 60.038044 22334.169 2013948.7
200000 300.63919 -9.0297919 550.45902 60.939218 60.038044 22943.626 2013948.7
Loop time of 81153.5 on 16 procs for 200000 steps with 149746 atoms
Pair time (%) = 52442.4 (64.6212)
Kspce time (%) = 23661.2 (29.1561)
Neigh time (%) = 1437.58 (1.77143)
Comm time (%) = 273.439 (0.33694)
Outpt time (%) = 9.62033 (0.0118545)
Other time (%) = 3329.26 (4.10243)
FFT time (% of Kspce) = 5540.4 (23.4156)
FFT Gflps 3d (1d only) = 6.54116 22.1496
Nlocal: 9359.12 ave 9468 max 9243 min
Histogram: 3 1 0 0 2 3 3 1 2 1
Nghost: 17619.2 ave 17753 max 17390 min
Histogram: 1 0 0 0 3 2 3 5 0 2
Neighs: 1.62895e+06 ave 1.6651e+06 max 1.59081e+06 min
Histogram: 1 1 1 1 3 3 2 2 1 1
FullNghs: 3.2579e+06 ave 3.32762e+06 max 3.18553e+06 min
Histogram: 1 2 0 2 3 2 1 3 1 1
Total # of neighbors = 52126428
Ave neighs/atom = 348.099
Neighbor list builds = 18643
Dangerous builds = 227
Thank you so much for you kind help!
Tao Wang / 王涛
Interdisciplinary Centre for Advanced Materials Simulation (ICAMS)
Ruhr-Universität Bochum
Stiepeler Str. 129
44801 Bochum
Germany
Tel: +49 234 32 29377
Fax: +49 234 32 14984