Dear Lammps developer,
I am working on reaxff, and using a model that the liquid film is confined between two solid substrates.
I plan to use fix qeq/reaxff to realize charge equilibrium in the system,
but I found the result was the liquid film became electronegativity, and the solid substrates became electropositivity. I wonder if the reason is I applied this command on the entire system. (the substrate is diamond, and there are many carbon atoms in the liquid film)
Because I only want to focus on the chemical reactions in liquid, the solids are just for providing the confinement and the atomic charge for solids is not so important. So I think if I can only apply fix qeq/reaxff on liquid.
I have tried it with 2 methods seperately, but quickly errors occured.
- Change the fix group-ID to liquid.
- In the required file, only list the “itype chi eta gamma” parameters for the atom types in liquid.
So I want to know if there are any methods to apply fix qeq/reaxff only on liquid in my system.
Best regards,
Xingyu Chen
First off, please always report which version of LAMMPS you are using and what platform you are running on.
what errors?
Method 1 should do what you ask for.
If all you care about is a confinement, you could also just use one of the “fix wall” commands and not have any atoms there at all.
1 Like
Dear Dr. Axel,
Thank you for your recommendation.
The version of LAMMPS is 29Sep2021,
and the platform is Linux with the CPU:
Intel Xeon Gold 6230, 20 core, 2.10 - 3.90 GHz Ă— 4
For method 1, the error is as follows:
---------------- Step 0 ----- CPU = 0.0000 (sec) ----------------
TotEng = 12866207.5605 KinEng = 13947.5083 Temp = 167.7879
PotEng = 12852260.0522 E_bond = 0.0000 E_angle = 0.0000
E_dihed = 0.0000 E_impro = 0.0000 E_vdwl = -4152597.8079
E_coul = 17004857.8600 E_long = 0.0000 Press = 354110.1749
ERROR on proc 0: Non-numeric atom coords - simulation unstable (…/domain_omp.cpp:58)
Last command: run 10000
ERROR on proc 1: Non-numeric atom coords - simulation unstable (…/domain_omp.cpp:58)
Last command: run 10000
…
(until proc 79)
By the way, the E_coul is so high and seems unreasonable. If I apply fix qeq/reaxff on the whole system, then the E_coul is -25738.7154.
Best regards,
Xingyu Chen
That would suggest that possibly your initial charge values for the boundary atoms is not suitable.
For a test, you should try running ReaxFF without fix qeq/reaxff (after disabling the check for it). It should give reasonable forces and energies.
Thanks for your suggestion. I tested it, and the E_coul became -6522.1657. I think it is more reasonable than if only using qeq on liquid.
I found an example in /example/reaxff/CHO. It set the charge for all atoms to 0.0 in the data file. So I think if I used fix qeq, then there is no problem even though I do not defined the initial charge.
As a test, I set the initial charge for all atoms in my model to 0.0, and only use fix qeq on liquid. But still I got the error as follow:
---------------- Step 0 ----- CPU = 0.0000 (sec) ----------------
TotEng = 15011070.5162 KinEng = 13947.5083 Temp = 167.7879
PotEng = 14997123.0079 E_bond = 0.0000 E_angle = 0.0000
E_dihed = 0.0000 E_impro = 0.0000 E_vdwl = -4152597.8078
E_coul = 19149720.8158 E_long = 0.0000 Press = 323664.4543
WARNING: Fix qeq/reaxff CG convergence failed after 10000 iterations at step 423 (…/fix_qeq_reaxff.cpp:746)
ERROR on proc 75: step 422: hbondchk failed: H=27 end(H)=746 str(H+1)=743
(…/reaxff_forces.cpp:119)
I just made a test setting up a layered system with QeQ only applied to the liquid part.
My conclusion of this has to be that your input or data file is not correct or perhaps your force field parameters are not suitable for your system or a combination of those.
Here is my example output with fix qeq/reaxff applied to only water atoms.
c_qC, c_qO, c_qH are the average charges per carbon, oxygen and hydrogen atom and q_qwater is the total(!) sum of the charge of the water atoms
Step E_vdwl E_coul TotEng Press c_qC c_qO c_qH c_qwater
0 -35338.558 -2759.8541 -37898.102 50416.859 0 -0.66641055 0.33320527 5.7176486e-15
100 -37883.181 -2859.0847 -38145.715 -26280.181 0 -0.69020417 0.34510208 -6.6613381e-16
200 -37846.078 -2895.9964 -38650.871 27161.248 0 -0.69905889 0.34952944 -3.6082248e-15
300 -37941.974 -2893.283 -39120.297 9712.0404 0 -0.698411 0.3492055 -2.1094237e-15
400 -38030.524 -2928.1035 -39516.664 -58165.64 0 -0.7067544 0.3533772 3.0531133e-15
500 -38056.588 -2953.9592 -39843.535 -31100.378 0 -0.71295992 0.35647996 -1.6653345e-16
600 -38238.921 -2953.3942 -40125.793 28090.752 0 -0.71282312 0.35641156 -1.9984014e-15
700 -38383.823 -2940.2657 -40355.241 -12575.521 0 -0.70967335 0.35483668 -4.8849813e-15
800 -38511.521 -2935.9466 -40561.871 -35956.335 0 -0.70863728 0.35431864 2.3314684e-15
900 -38607.152 -2916.8141 -40738.271 -34137.438 0 -0.70405351 0.35202676 -1.4432899e-15
1000 -38716.253 -2918.0475 -40890.174 -22292.369 0 -0.70435112 0.35217556 -1.9984014e-15
Here is the output for the same system but this time with QeQ applied to the entire system. As you can see there is some polarization of the carbon atoms, but it is two orders of magnitude lower than those of oxygen and hydrogen. And the total net charge water sub-system is rather low.
Step E_vdwl E_coul TotEng Press c_qC c_qO c_qH c_qwater
0 -35338.558 -2789.6115 -37927.859 50042.416 -0.0076939304 -0.66777634 0.34004331 0.92327165
100 -37880.2 -2868.242 -38164.322 -30498.784 -0.0048953914 -0.6887337 0.34828316 0.58744697
200 -37812.766 -2911.0883 -38659.648 25249.783 -0.0074116099 -0.69712942 0.354494 0.88939319
300 -37951.174 -2912.1628 -39126.097 11663.33 -0.0063776186 -0.69816284 0.35418352 0.76531423
400 -37973.655 -2928.9967 -39512.251 -30937.234 -0.0071390499 -0.7016157 0.35651909 0.85668599
500 -38102.645 -2966.0718 -39835.37 -54314.498 -0.0065098228 -0.71098353 0.36069962 0.78117874
600 -38203.929 -2968.284 -40104.685 7363.8564 -0.0062346795 -0.71172545 0.36085047 0.74816154
700 -38345.483 -2958.4214 -40334.553 -7316.7171 -0.0062631047 -0.70933748 0.35967922 0.75157257
800 -38420.205 -2947.1721 -40533.505 -12074.104 -0.0056098969 -0.70713272 0.35805428 0.67318763
900 -38549.409 -2936.4245 -40706.671 -31677.013 -0.0042349405 -0.70558448 0.35618019 0.50819286
1000 -38638.558 -2924.0753 -40856.447 -51806.948 -0.0043169277 -0.70256607 0.35473658 0.51803133
Loop time of 29.2468 on 1 procs for 1000 steps with 345 atoms
Dear Dr. Axel,
Thank you very much for your kind response and testing for me. Your answers really help me a lot. I think I have to carefully check my force field parameters and system.
Best regards,
Xingyu Chen