rohskopf / modecode

Massively parallel vibrational mode calculator.
21 stars 8 forks source link

Atomic displacement calculations are not generalized for all box types in LAMMPS. #1

Open rohskopf opened 3 years ago

rohskopf commented 3 years ago

In the LAMMPS compute and fix files, we use atomic displacements for some mode calculations, and this displacements are calculated via the following code snippet:

  // Compute atomic displacements.
  int iindx;
  for (int i=0; i<nlocal; i++){
    iindx = tag[i]-1;
    if (mask[i] & groupbit){ // This keeps the atoms in a particular group.
      for (int a=0; a<3; a++){
        u_p[i][a] = (x[i][a]-x0[tag[i]-1][a]);
        if (std::abs(u_p[i][a]) > domain->boxhi[a]/2.0 && u_p[i][a] > 0) u_p[i][a] -= domain->boxhi[a];
        else if (std::abs(u_p[i][a]) > domain->boxhi[a]/2.0 && u_p[i][a] < 0) u_p[i][a] += domain->boxhi[a];

      }
    }
  }

The way I've enforced periodic boundary conditions doesn't allow this to be applied to boxes whose lower bounds don't start at zero. For the systems I'm studying now, this is sufficient, but this should be fixed to generalize to all box shapes and sizes.