I want to compute the dipole moment of a spherical region in the bulk liquid water. I test two routes to obtain the dipole moment:
1 . Define a spherical region, and make a dynamic group to put all atoms in this spherical region volume into the group, then use “variable Mx atom x*q” and “compute reduce sum” to get dipole moment.
Use “chunk/atom” to define a spherical region, and compute the dipole moment by “compute dipole/chunk”.
My problem is those two routes can not give the same values of the dipole moment. In addition, as it is mentioned in the manual about “compute dipole/chunk” command, one may reset the image flags (e.g. to 0) before invoking this compute by using the “set image” command. However, there is no problem if I want to compute the total dipole moment of the whole bulk system, but it gives weird results for a spherical region.
you are likely not taking into account what it is exactly that compute dipole/chunk computes. specifically, there is:
For chunks with a net charge the resulting dipole is made position independent by subtracting the position vector of the center of mass or geometric center times the net charge from the computed dipole vector.
so your computation of: Sum_n q_ir_i is not what compute dipole/chunk is computing which is instead (by default): Sum_n q_ir_i - Sum_n q_i * (Sum_n m_ir_i)/Sum_n m_i
or Sum_n q_ir_i - Sum_n q_i * Sum_nr_i (when using the “geom” option).
when cutting a spherical region of water molecules, you will likely have a net charge. with water molecules, the cleaner solution in either case would be to make the group selection or chunk assignment such that whole molecules are selected where then the ambiguity and position dependence of Sum_n q_i*r_i due to a net charge is avoided by construction.
I still have one more question about the command “set image”. The image flag is reset before the “compute dipole/chunk” command, the value of the dipole moment is much greater than that without "set image"command. While there is no problem if I want to compute the total dipole moment of the whole system.
you must not blindly just reset all image flags to zero, especially not for molecules. in cases where a molecule stretches across a periodic boundary, you would break it into pieces as normally part of the molecule would have different image flags than other parts.
for neutral molecules in general, it would be best (for total dipole moment and per chunk/group/region dipole moment) to do the selection on just one reference atom of the molecule based on wrapped coordinates (i.e. without consideration of image flags) and then selecting the whole molecule as then the computation of the dipole moment should use unwrapped coordinates and will be position independent and more meaningful.
if there needs to be a reset of image flags, one might try writing out a data file and then processing the image flags in there with a custom script. or do it on a dump file and the first read the data file and then update the positions (and image flags) with read_data.
sorry, but what you are telling us is getting more inconsistent with every e-mail and thus leading nowhere.
you need to start over and provide a suitable, small but complete input deck and a specific explanation of where you see a problem and describe in detail how to reproduce it.
i don’t see any problems with my examples and cannot tell whether you are doing things correctly or if the situation is exactly as you claim it is.
Actually, the volume of the spherical region is smaller than the whole simulation box. The molecules in this sphere must be complete, none of them across the periodic boundary. So they all have the same image flags, “set image” should not affect the result.
And, the dipole moment should increase with increasing volume(M=P*V, M: dipole moment, P:polarization, V:volume). However, as I change the radius of the spherical region, the “compute dipole/moment” always give the same result, even the number after the decimal point is the same.
The image flag is reset before the “compute dipole/chunk” command, the value of the dipole moment is much greater than that without "set image"command. While there is no problem if I want to compute the total dipole moment of the whole system. The volume of the spherical region is smaller than the whole simulation box. The molecules in this sphere must be complete, none of them across the periodic boundary. So they all have the same image flags, “set image” should not affect the result.
The dipole moment should increase with increasing volume(M=P*V, M: dipole moment, P:polarization, V:volume). However, as I change the radius of the spherical region ${rc}, the “compute dipole/moment” always give the same result, even the number after the decimal point is the same.