Commit 10034ce3 authored by Axel Kohlmeyer's avatar Axel Kohlmeyer
Browse files

port support for scale[] factor with fix adapt to OPT and USER-OMP

parent 281ace32
Loading
Loading
Loading
Loading
+5 −2
Original line number Original line Diff line number Diff line
@@ -247,6 +247,7 @@ void PairEAMOpt::eval()
    if (EFLAG) {
    if (EFLAG) {
      double phi = ((coeff[3]*p + coeff[4])*p + coeff[5])*p + coeff[6];
      double phi = ((coeff[3]*p + coeff[4])*p + coeff[5])*p + coeff[6];
      if (rho[i] > rhomax) phi += fp[i] * (rho[i]-rhomax);
      if (rho[i] > rhomax) phi += fp[i] * (rho[i]-rhomax);
      phi *= scale[type[i]][type[i]];
      if (eflag_global) eng_vdwl += phi;
      if (eflag_global) eng_vdwl += phi;
      if (eflag_atom) eatom[i] += phi;
      if (eflag_atom) eatom[i] += phi;
    }
    }
@@ -273,6 +274,7 @@ void PairEAMOpt::eval()
    double tmpfz = 0.0;
    double tmpfz = 0.0;


    fast_gamma_t* _noalias tabssi = &tabss[itype1*ntypes*nr];
    fast_gamma_t* _noalias tabssi = &tabss[itype1*ntypes*nr];
    double* _noalias scale_i = scale[itype1+1]+1;


    for (jj = 0; jj < jnum; jj++) {
    for (jj = 0; jj < jnum; jj++) {
      j = jlist[jj];
      j = jlist[jj];
@@ -316,12 +318,13 @@ void PairEAMOpt::eval()
        // psip needs both fp[i] and fp[j] terms since r_ij appears in two
        // psip needs both fp[i] and fp[j] terms since r_ij appears in two
        //   terms of embed eng: Fi(sum rho_ij) and Fj(sum rho_ji)
        //   terms of embed eng: Fi(sum rho_ij) and Fj(sum rho_ji)
        //   hence embed' = Fi(sum rho_ij) rhojp + Fj(sum rho_ji) rhoip
        //   hence embed' = Fi(sum rho_ij) rhojp + Fj(sum rho_ji) rhoip
        // scale factor can be applied by thermodynamic integration


        double recip = 1.0/r;
        double recip = 1.0/r;
        double phi = z2*recip;
        double phi = z2*recip;
        double phip = z2p*recip - phi*recip;
        double phip = z2p*recip - phi*recip;
        double psip = fp[i]*rhojp + fp[j]*rhoip + phip;
        double psip = fp[i]*rhojp + fp[j]*rhoip + phip;
        double fpair = -psip*recip;
        double fpair = -scale_i[jtype]*psip*recip;


        tmpfx += delx*fpair;
        tmpfx += delx*fpair;
        tmpfy += dely*fpair;
        tmpfy += dely*fpair;
@@ -332,7 +335,7 @@ void PairEAMOpt::eval()
          ff[j].z -= delz*fpair;
          ff[j].z -= delz*fpair;
        }
        }


        if (EFLAG) evdwl = phi;
        if (EFLAG) evdwl = scale_i[jtype]*phi;


        if (EVFLAG) ev_tally(i,j,nlocal,NEWTON_PAIR,
        if (EVFLAG) ev_tally(i,j,nlocal,NEWTON_PAIR,
                             evdwl,0.0,fpair,delx,dely,delz);
                             evdwl,0.0,fpair,delx,dely,delz);
+4 −3
Original line number Original line Diff line number Diff line
@@ -203,7 +203,7 @@ void PairEAMOMP::eval(int iifrom, int iito, ThrData * const thr)
    if (EFLAG) {
    if (EFLAG) {
      phi = ((coeff[3]*p + coeff[4])*p + coeff[5])*p + coeff[6];
      phi = ((coeff[3]*p + coeff[4])*p + coeff[5])*p + coeff[6];
      if (rho[i] > rhomax) phi += fp[i] * (rho[i]-rhomax);
      if (rho[i] > rhomax) phi += fp[i] * (rho[i]-rhomax);
      e_tally_thr(this, i, i, nlocal, NEWTON_PAIR, phi, 0.0, thr);
      e_tally_thr(this, i, i, nlocal, NEWTON_PAIR, scale[type[i]][type[i]]*phi, 0.0, thr);
    }
    }
  }
  }


@@ -230,6 +230,7 @@ void PairEAMOMP::eval(int iifrom, int iito, ThrData * const thr)
    ztmp = x[i].z;
    ztmp = x[i].z;
    itype = type[i];
    itype = type[i];
    fxtmp = fytmp = fztmp = 0.0;
    fxtmp = fytmp = fztmp = 0.0;
    const double * _noalias const scale_i = scale[itype];


    jlist = firstneigh[i];
    jlist = firstneigh[i];
    jnum = numneigh[i];
    jnum = numneigh[i];
@@ -274,7 +275,7 @@ void PairEAMOMP::eval(int iifrom, int iito, ThrData * const thr)
        phi = z2*recip;
        phi = z2*recip;
        phip = z2p*recip - phi*recip;
        phip = z2p*recip - phi*recip;
        psip = fp[i]*rhojp + fp[j]*rhoip + phip;
        psip = fp[i]*rhojp + fp[j]*rhoip + phip;
        fpair = -psip*recip;
        fpair = -scale_i[jtype]*psip*recip;


        fxtmp += delx*fpair;
        fxtmp += delx*fpair;
        fytmp += dely*fpair;
        fytmp += dely*fpair;
@@ -285,7 +286,7 @@ void PairEAMOMP::eval(int iifrom, int iito, ThrData * const thr)
          f[j].z -= delz*fpair;
          f[j].z -= delz*fpair;
        }
        }


        if (EFLAG) evdwl = phi;
        if (EFLAG) evdwl = scale_i[jtype]*phi;
        if (EVFLAG) ev_tally_thr(this, i,j,nlocal,NEWTON_PAIR,
        if (EVFLAG) ev_tally_thr(this, i,j,nlocal,NEWTON_PAIR,
                                 evdwl,0.0,fpair,delx,dely,delz,thr);
                                 evdwl,0.0,fpair,delx,dely,delz,thr);
      }
      }