[lammps-users] Antw: Re: group/group for colloids

Dear Steve

Thank you very much.
Well, I am not sure whether I understood the first sentence correctly, but I clearly distinguish in the input between colloid1 and 2 as you can see below

#define groups

group liquid type 1
group colloid1 id 105134
group colloid2 id 105135
group colloid union colloid1 colloid2

corresponding to the log. output during run

105133 atoms in group liquid
1 atoms in group colloid1
1 atoms in group colloid2
2 atoms in group colloid

which looks fine. The problems are as mentioned in the last mail not the colloid/colloid interactions but the colloid/liquid ones.
I have done as you advised and printed out fpair and evwdl from compute() as well as fforce and phi from single() for the case SMALL_LARGE. The energies are in both cases the same, but the forces differ by a magnitude of two, as I have observed by comparison of the dump and compute_group_group output.

I doubt the forces from the group/group command deduced from single(), since they are unphysically large in x and y direction, where both particles lie on the z axis. Maybe there is a typo in the formula for SMALL-LARGE interaction, for instance fforce is multiplied with r fpair not. I have to check the derivative tomorrow, hopefully thats it. Do you have another idea?

Another thing I striked is in the routine pair.cpp at mix_energy() where for SIXTH
POWER mixing is written: (Version Okt21,2010)

value = 2.0 * sqrt(eps1*eps2) *pow(sig1,3.0) * pow(sig2,3.0) / (pow(sig1,6.0) * pow(sig2,6.0));

I think the formula should be the same as in the manual with + instead of * in the denominator

value = 2.0 * sqrt(eps1*eps2) *pow(sig1,3.0) * pow(sig2,3.0) / (pow(sig1,6.0) + pow(sig2,6.0));

Best regards

Steve Plimpton [email protected] 20.10.10 16.59 Uhr >>>
I don’t know. I think you are saying that the colloid1 and 2 groups
are just a single particle? I would debug this by using the
compute pair/local command and dumping out all pairwise forces
one by one, then summing them by hand. You could also put
print statements in the pair_colloid.cpp routine and see if the
answer in the compute() routine is the same as the one
in the single() routine. It sounds like you only have one
colloid-colloid interaction, so this should be easy to check.


I found an extra factor of "r" in the small-large formula in
the single() function of pair colloid which would affect the
compute group/group result (for force). This rings a faint
bell. I think it was a bug fix we made in the compute() routine
months ago, but must have neglected to make in the single().

Also, you're right about the sixthpower mixing formlula - that's
a bug also.

I'll post patches for both of these this AM. See if you get
good colloid agreement after that.