Commit a1f421cd authored by Aidan Thompson's avatar Aidan Thompson
Browse files

Moved compute_beta outside of main force loop

parent 6d84bd61
Loading
Loading
Loading
Loading
+25 −12
Original line number Diff line number Diff line
@@ -101,6 +101,7 @@ PairSNAP::PairSNAP(LAMMPS *lmp) : Pair(lmp)

  sna = NULL;

  beta_max = 0;
}

/* ---------------------------------------------------------------------- */
@@ -189,6 +190,15 @@ void PairSNAP::compute_regular(int eflag, int vflag)
  int newton_pair = force->newton_pair;
  class SNA* snaptr = sna[0];

  if (beta_max < list->inum) {
    memory->grow(beta,list->inum,ncoeff,"PairSNAP:beta");
    beta_max = list->inum;
  }

  // compute dE_i/dB_i = beta_i for all i in list

  compute_beta();

  numneigh = list->numneigh;
  firstneigh = list->firstneigh;

@@ -250,13 +260,9 @@ void PairSNAP::compute_regular(int eflag, int vflag)
    // compute Fij = dEi/dRj = -dEi/dRi 
    // add to Fi, subtract from Fj

    // compute dE_i/dB_i = beta_i

    compute_betai(ielem);

    // compute beta_i*Z_i = Y_i

    snaptr->compute_yi(beta);
    snaptr->compute_yi(beta[ii]);

    for (int jj = 0; jj < ninside; jj++) {
      int j = snaptr->inside[jj];
@@ -1294,15 +1300,23 @@ void PairSNAP::build_per_atom_arrays()
}

/* ----------------------------------------------------------------------
   compute beta_i 
   compute beta
------------------------------------------------------------------------- */

void PairSNAP::compute_betai(int ielem)
void PairSNAP::compute_beta()
{
  int i;
  int *type = atom->type;

  for (int ii = 0; ii < list->inum; ii++) {
    i = list->ilist[ii];
    const int itype = type[i];
    const int ielem = map[itype];
    double* coeffi = coeffelem[ielem];

    for (int k = 1; k <= ncoeff; k++)
    beta[k-1] = coeffi[k];
      beta[ii][k-1] = coeffi[k];
  }
}

/* ----------------------------------------------------------------------
@@ -1631,7 +1645,6 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename)
  memory->create(radelem,nelements,"pair:radelem");
  memory->create(wjelem,nelements,"pair:wjelem");
  memory->create(coeffelem,nelements,ncoeffall,"pair:coeffelem");
  memory->create(beta,ncoeffall,"pair:beta");

  // Loop over nelements blocks in the SNAP coefficient file

+3 −2
Original line number Diff line number Diff line
@@ -55,7 +55,7 @@ protected:
  void set_sna_to_shared(int snaid,int i);
  void build_per_atom_arrays();

  void compute_betai(int);
  void compute_beta();

  int schedule_user;
  double schedule_time_guided;
@@ -101,11 +101,12 @@ protected:
  double *radelem;              // element radii
  double *wjelem;               // elements weights
  double **coeffelem;           // element bispectrum coefficients
  double* beta;                 // beta for current atom
  double** beta;                // betas for all atoms in list
  int *map;                     // mapping from atom types to elements
  int twojmax, diagonalstyle, switchflag, bzeroflag;
  double rfac0, rmin0, wj1, wj2;
  int rcutfacflag, twojmaxflag; // flags for required parameters
  int beta_max;                 // length of beta
};

}