Hello Lammps Users,
Is it possible to achieve polarization of two specific parts of the simulation setup by ionic charges and then perform charge equilibration using fix qeq/reax? I have attached a schematic which shows what I want.
I tried to achieve this but invariably the charge that I put in gets overwritten to a negligible (close to zero) value right from the first step of charge equilibration (fix qeq/reax). Am I missing something here or is it something not possible with the current qeq/reax algorithm, unless we change something?
Thanks a lot,
Rensselaer Polytechnic Institute
In principle the QEq method doesn’t require charge neutrality, but in the implementation of the matrix inversion method (qeq/reax) it was found that non-charge neutrality leads to numerical instabilities causing charges to have difficulties to converge. What you can do is use one of the damped dynamics method (particularly qeq/dynamic) for charge equilibration. The downside is that you have prepare the “param.qeq” yourself, which shall be extracted from the ReaxFF forcefield file.
More to clarify. Qeq/reax neutralizes the whole system by setting all charges to zero… perhaps some modifications can be done to qeq/reax or qeq/shielded so that charge neutrality enforcement does not wipe the entire charges out.
Thanks a lot for your reply Ray. I was trying out fix qeq/dynamic, using a basic setup where I am assigning positive and negative charges to two blocks, so essentially the net charge is zero but there is a charge polarization. But using fix qeq/dynamic too the charge that I put in gets overwritten to a zero value right from the first step of charge equilibration.
I have attached the corresponding files with this email, in case you could give me a feedback if I am doing something wrong, or whether this is a limitation of fix qeq/dynamic too…
discharge.dis (1.75 KB)
ffield.reax (5.38 KB)
log.lammps (4 KB)
param.sili (75 Bytes)
In your input script, your are applying fix qeq to the whole simulation box, which is why your left/right-region charges get overwritten. Try applying the fix only to the central region.
In your input script, your are applying fix qeq to the whole simulation
box, which is why your left/right-region charges get overwritten. Try
applying the fix only to the central region.
good idea, but i think while this will keep those charges untouched, it
will also not include them in the qeq calculation, which is not desired.
from a quick glance at the source code it looks that to have charges
considered, but not changed during qeq, will likely require some
programming, where a second group is introduced that flags the "do not
change the charge of these atoms" condition.
Arthur and Axel provided good suggestions.
The reason why your charged blocks became charge neutral after one charge equilibration step wasn’t because of fix qeq/dynamic wiping out the charges, but simply due to charge equilibration. Please remember that QEq works by equalizing the electronegativity (or derivative of energy wrt charge) on atoms. If your system has atoms with the same electronegativities, no matter what their initial charges are, charges are re-distributed so that partial charges on atoms are same throughout the system. This was what happened to your system. Since charges are balanced, re-distribution of charges meant neutral charge throughout the system. If you have positive net charge, this net charge will be equally divided to all the atoms. This is what happens when your system has atoms with the same electronegativity.
Now, if you create “electrodes” with different electronegativities, e.g. one block with Si atoms and the other with O atoms, you should see positive charges automatically created for the Si block and negative charges for the O block. If you further exclude these two blocks from time-integration, you effectively create a polarization throughout the system. That may be a word-around for you.
Lastly, the purpose of polarization throughout the system is to create an electric field. If this is indeed your purpose, it would be much simpler to use fix efield to create the external field.
Thank you Arthur, Axel and Ray for your wonderful suggestions…