[lammps-users] fix rigid temperature compute

Hello,

I believe that for fix rigid, the DOF computation in two dimensions is incorrect. For each rigid body, 2N-3 degrees should be subtracted.

something like this in the code:

   int n = 0;
   if (domain->dimension == 3) {
     for (int ibody = 0; ibody < nbody; ibody++) {
       if (nall[ibody] > 2) n += 3*nall[ibody] - 6;
       else if (nall[ibody] == 2) n++;
     }
   } else if (domain->dimension == 2) {
     for (int ibody = 0; ibody < nbody; ibody++) n += 2*nall[ibody] - 3;
   }

Alternatively, this can be corrected by a compute_modify.

Yours,
Tony

Hi Tony. Yes, I think you are correct. Another way to code it would be something like this:

  // remove d*N - d(d+1)/2 dof for each rigid body if more than 2 atoms in igroup
  // remove an additional dof for each diatomic rigid body in igroup for 3D simulations
  int n = 0;
  for (int ibody = 0; ibody < nbody; ibody++) {
    if (nall[ibody] > 1) n += domain->dimension*nall[ibody] - domain->dimension*(domain->dimension + 1)/2;
    if ((nall[ibody] == 2) and (domain->dimension == 3)) n++;
  }

(See discussion here:

http://en.wikipedia.org/wiki/Degrees_of_freedom_(engineering) ).

Your version needs to account for the (somewhat trivial) case of 1 atom in the rigid body in 2D simulations, where the dof should remain unchanged.

Paul

Posted the patch for this.

Thanks,
Steve