ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

move eflag

Browse Source
This commit is contained in:
Eddy Barraud
2024-06-06 10:10:36 +02:00
parent 8eefc0d305
commit c1db6a50ec
1 changed files with 43 additions and 47 deletions
Download Patch File Download Diff File Expand all files Collapse all files

76
lib/gpu/lal_dpd_charged.cu
View File

@ -185,7 +185,7 @@ __kernel void k_dpd_charged(const __global numtyp4 *restrict x_,
atom_info(t_per_atom,ii,tid,offset);
__local numtyp sp_cl[4];
local_allocate_store_charge();
///local_allocate_store_charge();
sp_cl[0]=sp_cl_in[0];
sp_cl[1]=sp_cl_in[1];
@ -200,7 +200,7 @@ __kernel void k_dpd_charged(const __global numtyp4 *restrict x_,
acctyp e_coul, energy, virial[6];
if (EVFLAG) {
energy=(acctyp)0;
e_coul=(acctyp)0
e_coul=(acctyp)0;
for (int i=0; i<6; i++) virial[i]=(acctyp)0;
}
@ -270,14 +270,22 @@ __kernel void k_dpd_charged(const __global numtyp4 *restrict x_,
// drag force = -gamma * wd^2 * (delx dot delv) / r
// random force = sigma * wd * rnd * dtinvsqrt;
if (!tstat_only) force_dpd = coeff[mtype].x*wd;
force_dpd -= coeff[mtype].y*wd*wd*dot*rinv;
force_dpd *= factor_dpd;
force_dpd += factor_sqrt*coeff[mtype].z*wd*randnum*dtinvsqrt;
force_dpd *=rinv;
if (EVFLAG && eflag) {
// unshifted eng of conservative term:
// evdwl = -a0[itype][jtype]*r * (1.0-0.5*r/cut[itype][jtype]);
// eng shifted to 0.0 at cutoff
numtyp e = (numtyp)0.5*coeff[mtype].x*coeff[mtype].w * wd*wd;
energy += factor_dpd*e;
}
}// if cut_dpdsq
// apply Slater electrostatic force if distance below Slater cutoff
// and the two species have a slater coeff
// cutsq[mtype].w -> Coulombic squared cutoff
@ -298,29 +306,20 @@ __kernel void k_dpd_charged(const __global numtyp4 *restrict x_,
force_coul = prefactor*(_erfc + EWALD_F*grij*expm2-slater_term);
if (factor_coul > (numtyp)0) force_coul -= factor_coul*prefactor*((numtyp)1.0-slater_term);
force_coul *= r2inv;
if (EVFLAG && eflag) {
numtyp e_slater = ((numtyp)1.0 + rlamdainv)*exprlmdainv;
numtyp e = prefactor*(_erfc-e_slater);
if (factor_coul > (numtyp)0) e -= factor_coul*prefactor*((numtyp)1.0 - e_slater);
e_coul += e;
}
} // if cut_coulsq
numtyp force = force_coul + force_dpd;
f.x += delx*force;
f.y += dely*force;
f.z += delz*force;
if (EVFLAG && eflag) {
if (rsq < cutsq[mtype].y && r < EPSILON) {
// unshifted eng of conservative term:
// evdwl = -a0[itype][jtype]*r * (1.0-0.5*r/cut[itype][jtype]);
// eng shifted to 0.0 at cutoff
numtyp e = (numtyp)0.5*coeff[mtype].x*coeff[mtype].w * wd*wd;
energy+=factor_dpd*e;
}
if ( cutsq[mtype].w != 0.0 && rsq < cutsq[mtype].w){
numtyp e_slater = ((numtyp)1.0 + rlamdainv)*exprlmdainv;
numtyp e = prefactor*(_erfc-e_slater);
if (factor_coul > (numtyp)0) e -= factor_coul*prefactor*((numtyp)1.0 - e_slater);
e_coul += e;
}
}
if (EVFLAG && vflag) {
virial[0] += delx*delx*force;
virial[1] += dely*dely*force;
@ -330,8 +329,7 @@ __kernel void k_dpd_charged(const __global numtyp4 *restrict x_,
virial[5] += dely*delz*force;
}
}
} // if cutsq
} // for nbor
} // if ii
@ -459,14 +457,22 @@ __kernel void k_dpd_charged_fast(const __global numtyp4 *restrict x_,
// drag force = -gamma * wd^2 * (delx dot delv) / r
// random force = sigma * wd * rnd * dtinvsqrt;
if (!tstat_only) force_dpd = coeff[mtype].x*wd;
force_dpd -= coeff[mtype].y*wd*wd*dot*rinv;
force_dpd *= factor_dpd;
force_dpd += factor_sqrt*coeff[mtype].z*wd*randnum*dtinvsqrt;
force_dpd *=rinv;
if (EVFLAG && eflag) {
// unshifted eng of conservative term:
// evdwl = -a0[itype][jtype]*r * (1.0-0.5*r/cut[itype][jtype]);
// eng shifted to 0.0 at cutoff
numtyp e = (numtyp)0.5*coeff[mtype].x*coeff[mtype].w * wd*wd;
energy += factor_dpd*e;
}
}// if cut_dpdsq
// apply Slater electrostatic force if distance below Slater cutoff
// and the two species have a slater coeff
// cutsq[mtype].w -> Coulombic squared cutoff
@ -487,29 +493,20 @@ __kernel void k_dpd_charged_fast(const __global numtyp4 *restrict x_,
force_coul = prefactor*(_erfc + EWALD_F*grij*expm2-slater_term);
if (factor_coul > (numtyp)0) force_coul -= factor_coul*prefactor*((numtyp)1.0-slater_term);
force_coul *= r2inv;
if (EVFLAG && eflag) {
numtyp e_slater = ((numtyp)1.0 + rlamdainv)*exprlmdainv;
numtyp e = prefactor*(_erfc-e_slater);
if (factor_coul > (numtyp)0) e -= factor_coul*prefactor*((numtyp)1.0 - e_slater);
e_coul += e;
}
} // if cut_coulsq
numtyp force = force_coul + force_dpd;
f.x += delx*force;
f.y += dely*force;
f.z += delz*force;
if (EVFLAG && eflag) {
if (rsq < cutsq[mtype].y && r < EPSILON) {
// unshifted eng of conservative term:
// evdwl = -a0[itype][jtype]*r * (1.0-0.5*r/cut[itype][jtype]);
// eng shifted to 0.0 at cutoff
numtyp e = (numtyp)0.5*coeff[mtype].x*coeff[mtype].w * wd*wd;
energy+=factor_dpd*e;
}
if ( cutsq[mtype].w != 0.0 && rsq < cutsq[mtype].w){
numtyp e_slater = ((numtyp)1.0 + rlamdainv)*exprlmdainv;
numtyp e = prefactor*(_erfc-e_slater);
if (factor_coul > (numtyp)0) e -= factor_coul*prefactor*((numtyp)1.0 - e_slater);
e_coul += e;
}
}
if (EVFLAG && vflag) {
virial[0] += delx*delx*force;
virial[1] += dely*dely*force;
@ -519,8 +516,7 @@ __kernel void k_dpd_charged_fast(const __global numtyp4 *restrict x_,
virial[5] += dely*delz*force;
}
}
} // if cutsq
} // for nbor
} // if ii