[gmx-developers] virial symmetry

[Berk Hess]

Hi,

I fixed it.
Note that this bug has no effect on the results (except for those 
off-diagonal virial components),
since using fully anisotropic pressure coupling with walls anyhow makes 
no sense.

Berk

Berk Hess wrote:
> Ah, thanks!
>
> That is indeed incorrect.
> The same correction should also be applied for x and y on the lower wall.
>
> But it is not really Newton's 3rd law here, but the Gromacs single sum 
> virial trick.
>
> PS in your last line [XX][XX] should be [XX][ZZ]
> (but I should check the order of the two components in the virial).
>
> Berk
>
> Alexander Herz wrote:
>> Hi,
>>
>> no I haven't changed the walls.
>>
>> the virial contains : r_x*F_z, r_y*F_z and r_z*F_z in components other
>> then [ZZ][ZZ]
>> But the only correction made to the virial for the walls is for the
>> zz,zz component (last line of the wall code you copied)!
>> So (unlike for the rest of the pair potentials like lennard jones,
>> estat, etc.) the wall contribution to the virial via the atom forces is
>> position dependent!
>>
>> To reproduce the problem:
>> make a small sim with 2 walls which are identical except density (make
>> em differ by a factor of 10 or so)
>> and then look at the average virial at the end of the md.log, it will be
>> asymmetric (so theta_ij!=theta_ji).
>>
>> to fix the problem I'd suggest:
>>
>> inside the force calculating loop:
>> instead of
>>
>> Fwall[w]+=(-F);
>>
>> do:
>> //force MUST act on the same x/y coord as the particle otherwise newtons
>> 3rd law is violated
>> fr->vir_wall[XX]+=-0.5*x[i][XX]*(-F);
>> fr->vir_wall[YY]+=-0.5*x[i][YY]*(-F);
>> fr->vir_wall[XX]+=-0.5*box_zz*(-F);
>>
>> and inside sim_util.c:
>>
>> vir_part[XX][XX]+=fr->vir_wall[XX];
>> vir_part[YY][ZZ]+=fr->vir_wall[YY];
>> vir_part[ZZ][ZZ]+=fr->vir_wall[ZZ];
>>
>> Alex
>>
>>
>> Berk Hess wrote:
>>  
>>> Hi,
>>>
>>> At the end of do_walls it says:
>>>
>>>  /* Since the lower wall is at z=0 it does not contribute
>>>   * to the single sum virial.
>>>   */
>>>  if (nwall == 2)
>>>    fr->vir_wall_zz = -0.5*box_zz*Fwall[1];
>>>
>>> This line should ensure that the virial only depends on the distance
>>> between the particles and the upper wall.
>>>
>>> I don't understand what you mean with asymmetric.
>>> The wall only acts on the z-component of the forces and both the force
>>> and distance have no x and y
>>> components, therefore the wall forces only occur in the zz element of
>>> the virial tensor.
>>>
>>> Have you modified the code such that this correction no longer holds?
>>>
>>> Berk
>>>
>>> Alexander Herz wrote:
>>>    
>>>> Having looked a bit closer into this, it seams that the definition of
>>>> the virial used in gromacs
>>>> (using absolute particle coords) is obviously only translation 
>>>> invariant
>>>> as long newtons 3rd law
>>>> is obeyed inside the system and no external forces act on the 
>>>> system (so
>>>> all forces sum to zero).
>>>> Apparently this condition is violated when using 2 walls with 
>>>> different
>>>> parameters like density???
>>>>
>>>> Alex
>>>>
>>>> Alexander Herz wrote:
>>>>  
>>>>      
>>>>> Hi guys,
>>>>>
>>>>> I'm having some problems with the walls feature (new in gromacs 4.0).
>>>>>
>>>>> If I understand correctly then I do get a symmetric virial tensor
>>>>> even for
>>>>> systems which are inhomogenous as long only pair potentials are used
>>>>>
>>>>> (Theta=-0.5*sum(r_ij*F_ij) where r_ij is r_j-r_i (so the distance
>>>>> between to particles acting on each other)
>>>>> and F_ij is the force acting on them;with F_ij=-F_ji the 
>>>>> contributions
>>>>> cancel out and the viral is symmetric )
>>>>>
>>>>> Now if I use 2 walls with different density I get an asymmetric 
>>>>> virial
>>>>> and looking at the wall's source
>>>>> it seams that for walls not r_ij is used but absolute wall coordinate
>>>>> (box_zz for wall2 and 0 for wall 1
>>>>> giving a zero contribution to the viral for wall1). Also for the
>>>>> individual atoms it seams the absolute position
>>>>> is used rather than r_ij (looking at calc_vir() in calcvir.c).
>>>>>
>>>>> So why do the walls use the absolute z coordinate instead of
>>>>> particle_z-wall_z and why
>>>>> do I get an asymmetric virial for walls with different density but 
>>>>> not
>>>>> for the same
>>>>> system without walls which is still asymmetric in z (some crystal
>>>>> occupies the lower z half off the system, the upper
>>>>> half is filled with water).
>>>>>
>>>>> Thx,
>>>>> Alex
>>>>>
>>>>> _______________________________________________
>>>>> gmx-developers mailing list
>>>>> gmx-developers at gromacs.org
>>>>> http://www.gromacs.org/mailman/listinfo/gmx-developers
>>>>> Please don't post (un)subscribe requests to the list. Use the www
>>>>> interface or send it to gmx-developers-request at gromacs.org.
>>>>>               
>>>> _______________________________________________
>>>> gmx-developers mailing list
>>>> gmx-developers at gromacs.org
>>>> http://www.gromacs.org/mailman/listinfo/gmx-developers
>>>> Please don't post (un)subscribe requests to the list. Use the www
>>>> interface or send it to gmx-developers-request at gromacs.org.
>>>>         
>>> _______________________________________________
>>> gmx-developers mailing list
>>> gmx-developers at gromacs.org
>>> http://www.gromacs.org/mailman/listinfo/gmx-developers
>>> Please don't post (un)subscribe requests to the list. Use the www
>>> interface or send it to gmx-developers-request at gromacs.org.
>>>     
>>
>> _______________________________________________
>> gmx-developers mailing list
>> gmx-developers at gromacs.org
>> http://www.gromacs.org/mailman/listinfo/gmx-developers
>> Please don't post (un)subscribe requests to the list. Use the www 
>> interface or send it to gmx-developers-request at gromacs.org.
>>   
>
> _______________________________________________
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the www 
> interface or send it to gmx-developers-request at gromacs.org.