Unverified Commit 245f949e authored by Axel Kohlmeyer's avatar Axel Kohlmeyer
Browse files

functions theta(), theta_mod(), and cutoff() are never used in pair style bop

parent e0ee0986
Loading
Loading
Loading
Loading
+0 −159
Original line number Diff line number Diff line
@@ -1035,145 +1035,6 @@ void PairBOP::gneigh()

/* ---------------------------------------------------------------------- */

void PairBOP::theta()
{
  int i,j,ii,jj,kk;
  int itype,jtype,i12;
  int temp_ij,temp_ik,temp_ijk;
  int n,nlocal,nall,ks;
  int nlisti;
  int *ilist;
  int *iilist;
  int **firstneigh;
  double rj2,rk2,rsq,ps;
  double rj1k1,rj2k2;
  double **x = atom->x;
  int *type = atom->type;

  nlocal = atom->nlocal;
  nall = nlocal+atom->nghost;
  ilist = list->ilist;
  firstneigh = list->firstneigh;
  if(update_list!=0)
    memory_theta_grow();
  else
    memory_theta_create();
  for (ii = 0; ii < nall; ii++) {
    if(ii<nlocal)
      i= ilist[ii];
    else
      i=ii;
    itype = map[type[i]]+1;

    iilist=firstneigh[i];
    nlisti=BOP_total[i];
    for(jj=0;jj<nlisti;jj++) {
      temp_ij=BOP_index[i]+jj;
      j=iilist[neigh_index[temp_ij]];
      jtype = map[type[j]]+1;
      if(itype==jtype)
        i12=itype-1;
      else if(itype<jtype)
        i12=itype*bop_types-itype*(itype+1)/2+jtype-1;
      else
        i12=jtype*bop_types-jtype*(jtype+1)/2+itype-1;
      if(i12>=npairs) {
        error->one(FLERR,"Too many atom pairs for pair bop");
      }
      disij[0][temp_ij]=x[j][0]-x[i][0];
      disij[1][temp_ij]=x[j][1]-x[i][1];
      disij[2][temp_ij]=x[j][2]-x[i][2];
      rsq=disij[0][temp_ij]*disij[0][temp_ij]
          +disij[1][temp_ij]*disij[1][temp_ij]
          +disij[2][temp_ij]*disij[2][temp_ij];
      rij[temp_ij]=sqrt(rsq);
      if(rij[temp_ij]<=rcut[i12])
        neigh_flag[temp_ij]=1;
      else
        neigh_flag[temp_ij]=0;
      if(rij[temp_ij]<=rcut3[i12])
        neigh_flag3[temp_ij]=1;
      else
        neigh_flag3[temp_ij]=0;
      ps=rij[temp_ij]*rdr[i12]+1.0;
      ks=(int)ps;
      if(nr-1<ks)
        ks=nr-1;
      ps=ps-ks;
      if(ps>1.0)
        ps=1.0;
      betaS[temp_ij]=((pBetaS3[i12][ks-1]*ps+pBetaS2[i12][ks-1])*ps+pBetaS1[i12][ks-1])*ps+pBetaS[i12][ks-1];
      dBetaS[temp_ij]=(pBetaS6[i12][ks-1]*ps+pBetaS5[i12][ks-1])*ps
          +pBetaS4[i12][ks-1];
      betaP[temp_ij]=((pBetaP3[i12][ks-1]*ps+pBetaP2[i12][ks-1])*ps
          +pBetaP1[i12][ks-1])*ps+pBetaP[i12][ks-1];
      dBetaP[temp_ij]=(pBetaP6[i12][ks-1]*ps+pBetaP5[i12][ks-1])*ps
          +pBetaP4[i12][ks-1];
      repul[temp_ij]=((pRepul3[i12][ks-1]*ps+pRepul2[i12][ks-1])*ps
          +pRepul1[i12][ks-1])*ps+pRepul[i12][ks-1];
      dRepul[temp_ij]=(pRepul6[i12][ks-1]*ps+pRepul5[i12][ks-1])*ps
          +pRepul4[i12][ks-1];
    }
  }
  for (ii = 0; ii < nall; ii++) {
    n=0;
    if(ii<nlocal)
      i= ilist[ii];
    else
      i=ii;
    iilist=firstneigh[i];
    nlisti=BOP_total[i];
    for(jj=0;jj<nlisti;jj++) {
      temp_ij=BOP_index[i]+jj;
      j=iilist[neigh_index[temp_ij]];
      rj2=rij[temp_ij]*rij[temp_ij];
      for(kk=jj+1;kk<nlisti;kk++) {
        if(cos_index[i]+n>=cos_total) {
          error->one(FLERR,"Too many atom triplets for pair bop");
        }
        temp_ik=BOP_index[i]+kk;
        temp_ijk=cos_index[i]+n;
        if(temp_ijk>=cos_total) {
          error->one(FLERR,"Too many atom triplets for pair bop");
        }
        rk2=rij[temp_ik]*rij[temp_ik];
        rj1k1=rij[temp_ij]*rij[temp_ik];
        rj2k2=rj1k1*rj1k1;
        if(temp_ijk>=cos_total) {
          error->one(FLERR,"Too many atom triplets for pair bop");
        }
        cosAng[temp_ijk]=(disij[0][temp_ij]*disij[0][temp_ik]+disij[1][temp_ij]
            *disij[1][temp_ik]+disij[2][temp_ij]*disij[2][temp_ik])/rj1k1;
        dcAng[temp_ijk][0][0]=(disij[0][temp_ik]*rj1k1-cosAng[temp_ijk]
              *disij[0][temp_ij]*rk2)/(rj2k2);
        dcAng[temp_ijk][1][0]=(disij[1][temp_ik]*rj1k1-cosAng[temp_ijk]
            *disij[1][temp_ij]*rk2)/(rj2k2);
        dcAng[temp_ijk][2][0]=(disij[2][temp_ik]*rj1k1-cosAng[temp_ijk]
            *disij[2][temp_ij]*rk2)/(rj2k2);
        dcAng[temp_ijk][0][1]=(disij[0][temp_ij]*rj1k1-cosAng[temp_ijk]
            *disij[0][temp_ik]*rj2)/(rj2k2);
        dcAng[temp_ijk][1][1]=(disij[1][temp_ij]*rj1k1-cosAng[temp_ijk]
            *disij[1][temp_ik]*rj2)/(rj2k2);
        dcAng[temp_ijk][2][1]=(disij[2][temp_ij]*rj1k1-cosAng[temp_ijk]
            *disij[2][temp_ik]*rj2)/(rj2k2);
        n++;
      }
    }
  }
}

/* ---------------------------------------------------------------------- */

void PairBOP::theta_mod()
{
  if(update_list!=0)
    memory_theta_grow();
  else
    memory_theta_create();
}

/* ---------------------------------------------------------------------- */

/*  The formulation differs slightly to avoid negative square roots
    in the calculation of Sigma^(1/2) of (a) Eq. 6 and (b) Eq. 11 */

@@ -5556,26 +5417,6 @@ void _noopt PairBOP::read_table(char *filename)
  }
}

/* ---------------------------------------------------------------------- */

double PairBOP::cutoff(double rp,double vrcut,int mode,double r)
{
  double tmp,tmp_beta,tmp_alpha,cut_store;

  if(mode==1) {
    tmp=(rsmall-rbig)*(r-rp)/(vrcut-rp)+rbig;
    cut_store=(erfc(tmp)-erfc(rsmall))/(erfc(rbig)-erfc(rsmall));
  }
  else {
    tmp_beta=log(log(rbig)/log(rsmall))/log(rp/vrcut);
    tmp_alpha=-log(rbig)/pow(rp,tmp_beta);
    cut_store=(exp(-tmp_alpha*pow(r,tmp_beta))-exp(-tmp_alpha*pow(vrcut
        ,tmp_beta)))/(exp(-tmp_alpha*pow(rp,tmp_beta))-exp(-tmp_alpha
        *pow(vrcut,tmp_beta)));
  }
  return(cut_store);
}

/* ----------------------------------------------------------------------
   memory usage of local atom-based arrays
------------------------------------------------------------------------- */
+0 −3
Original line number Diff line number Diff line
@@ -199,14 +199,11 @@ class PairBOP : public Pair {
  void setPrepul();
  void setSign();
  void gneigh();
  void theta();
  void theta_mod();
  double sigmaBo(int, int);
  double PiBo(int, int);
  void memory_theta_create();
  void memory_theta_destroy();
  void memory_theta_grow();
  double cutoff(double, double, int, double);

  void read_table(char *);
  void allocate();