[lammps-users] Bad principal moments

Hello everyone,

I am trying to use lammps to solve problems in contact mechanics. My system consists hemispherical asperity placed on top of a flat plate and this asperity will be deformed by a rigid plane. But first I want to apply the following heat treatment to this system: The system is to be maintained at constant temperature of 1350K for 0.25ns, then temperature drops to 1000K and is held constante for another 0.25ns; temperature drops again to 700K and is held constant for 0.25ns; new drop to 300K and another 0.25ns.

With the input file (file attached) I am using the simulation starts OK at 1350K but when the temperature drops to 1000K I get the error:

" Bad principal moments"

Can someone please help?



in.anneal (1.6 KB)

The message is coming from fix rigid, and has to do with
the computation being done to keep the plane of atoms rigid.

Can you see if you can reproduce this with as small and short
simulation? Perhaps the error will always happen when you
make the transition to the 2nd run ...


Tried that. The error appears whenever at the first step of the 2nd run. I found a solution for this problem but I don't know if it is acceptable:

When I use fix setforce 0.0 0.0 0.0 on the rigid plane and exclude this group from the fix NVE during the annealing everything is ok.

Now I've come across another problem when the deformation part of the simulation begins (when I start to move the rigid plane towards the asperity):

If I use the fix rigid on the rigid plane during the deformation, the plane can't deform the asperity and starts to rotate about its center of mass. To move the rigid plane I've used the command:

displace_atoms rigid_plane move 0 -0.1 0 units box

inside a loop

Solution I've found:

I don't include the rigid plane atoms in the group of atoms controled by fix NVE command. Also I zero the angular momentum of this group of atoms and make sure that force acting on them can only act on the movement direction (negative y direction)

fix HELP rigid_plane momentum 1 angular
fix LINE rigid_plane lineforce 0.0 -1.0 0.0

Does this work-arrounds make any sense? Is there any physical basis for this?

Thanks for the help


Steve Plimpton wrote:

Yes, you don't want to have fix nve applied to the rigid plane. This is
b/c fix rigid does time integration itself, so you don't want a 2nd integration
applied. However, it sounds like what you really want is to just
move the plane of atoms into the asperity. In that case, fix rigid isn't
what you want, since it will allow the rigid body to rotate. Rather you
want to just set the initial V of the plane downward and have no fix nve
on those atoms. It will move steadily down forever.

Or if you want the asperity to push back and slow down
the plane, then use fix aveforce on the plane with fix nve. That will keep
the plane together, but push it back in Z.


Yes I just want to move the plane into the asperity at a fixed velocity, just like in a compression test. But instead of fixing the velocity of the rigid plane I am using:

variable i loop 160
label loop
thermo 100
displace_atoms rigid_plane move 0 -0.1 0 units box
run 100 post no
next i
jump in.indent loop

So I have at least 100 timesteps for the system to equilibrate at each displacement of th rigid plane. In the end that I wanto to measure the temperature evolution of the asperity, the forces of acting on the rigid plane, the stresses at contact and the contact area.

So is it ok to:

fix HELP rigid_plane momentum 1 angular
fix LINE rigid_plane lineforce 0.0 -1.0 0.0

on the rigid plane and not include them in the fix nve?

Thanks again


Steve Plimpton wrote: