From 9ffc09f497322c09ef531f1eae7fcf6393a36052 Mon Sep 17 00:00:00 2001 From: Rafael Ravedutti Date: Thu, 7 Jul 2022 00:47:38 +0200 Subject: [PATCH] Add DEM kernel to parameter options Signed-off-by: Rafael Ravedutti --- lammps/force_dem.c | 4 ++-- lammps/includes/util.h | 1 + lammps/main.c | 53 +++++++++++++++++++++--------------------- lammps/util.c | 2 ++ 4 files changed, 31 insertions(+), 29 deletions(-) diff --git a/lammps/force_dem.c b/lammps/force_dem.c index 6ba1c70..698c587 100644 --- a/lammps/force_dem.c +++ b/lammps/force_dem.c @@ -32,7 +32,7 @@ // TODO: Joint common files for gromacs and lammps variants #include "../gromacs/includes/simd.h" -double computeForceDEMFullNeigh_plain_c(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *stats) { +double computeForceDemFullNeigh(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *stats) { int Nlocal = atom->Nlocal; int* neighs; MD_FLOAT k_s = param->k_s; @@ -153,7 +153,7 @@ double computeForceDEMFullNeigh_plain_c(Parameter *param, Atom *atom, Neighbor * return E-S; } -double computeForceDEMHalfNeigh(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *stats) { +double computeForceDemHalfNeigh(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *stats) { int Nlocal = atom->Nlocal; int* neighs; #ifndef EXPLICIT_TYPES diff --git a/lammps/includes/util.h b/lammps/includes/util.h index 5adabe2..e317afb 100644 --- a/lammps/includes/util.h +++ b/lammps/includes/util.h @@ -44,6 +44,7 @@ #define FF_LJ 0 #define FF_EAM 1 +#define FF_DEM 2 extern double myrandom(int*); extern void random_reset(int *seed, int ibase, double *coord); diff --git a/lammps/main.c b/lammps/main.c index 3551ce1..03760bf 100644 --- a/lammps/main.c +++ b/lammps/main.c @@ -49,6 +49,8 @@ extern double computeForceLJFullNeigh_plain_c(Parameter*, Atom*, Neighbor*, Stat extern double computeForceLJFullNeigh_simd(Parameter*, Atom*, Neighbor*, Stats*); extern double computeForceLJHalfNeigh(Parameter*, Atom*, Neighbor*, Stats*); extern double computeForceEam(Eam*, Parameter*, Atom*, Neighbor*, Stats*); +extern double computeForceDemFullNeigh(Parameter*, Atom*, Neighbor*, Stats*); +extern double computeForceDemHalfNeigh(Parameter*, Atom*, Neighbor*, Stats*); #ifdef USE_SIMD_KERNEL # define KERNEL_NAME "SIMD" @@ -121,14 +123,29 @@ void finalIntegrate(Parameter *param, Atom *atom) { } void printAtomState(Atom *atom) { - printf("Atom counts: Natoms=%d Nlocal=%d Nghost=%d Nmax=%d\n", - atom->Natoms, atom->Nlocal, atom->Nghost, atom->Nmax); + printf("Atom counts: Natoms=%d Nlocal=%d Nghost=%d Nmax=%d\n", atom->Natoms, atom->Nlocal, atom->Nghost, atom->Nmax); + // int nall = atom->Nlocal + atom->Nghost; + // for (int i=0; ix[i], atom->y[i], atom->z[i]); + // } +} - /* int nall = atom->Nlocal + atom->Nghost; */ +double computeForce(Eam *eam, Parameter *param, Atom *atom, Neighbor *neighbor, Stats *stats) { + if(param->force_field == FF_EAM) { + return computeForceEam(eam, param, atom, neighbor, stats); + } else if(param->force_field == FF_DEM) { + if(param->half_neigh) { + return computeForceDemHalfNeigh(param, atom, neighbor, stats); + } else { + return computeForceDemFullNeigh(param, atom, neighbor, stats); + } + } - /* for (int i=0; ix[i], atom->y[i], atom->z[i]); */ - /* } */ + if(param->half_neigh) { + return computeForceLJHalfNeigh(param, atom, neighbor, stats); + } + + return computeForceLJFullNeigh(param, atom, neighbor, stats); } int main(int argc, char** argv) { @@ -208,7 +225,7 @@ int main(int argc, char** argv) { printf("MD Bench: A minimalistic re-implementation of miniMD\n"); printf(HLINE); printf("-p : file to read parameters from (can be specified more than once)\n"); - printf("-f : force field (lj or eam), default lj\n"); + printf("-f : force field (lj, eam or dem), default lj\n"); printf("-i : input file with atom positions (dump)\n"); printf("-e : input file for EAM\n"); printf("-n / --nsteps : set number of timesteps for simulation\n"); @@ -231,16 +248,8 @@ int main(int argc, char** argv) { #if defined(MEM_TRACER) || defined(INDEX_TRACER) traceAddresses(¶m, &atom, &neighbor, n + 1); #endif - if(param.force_field == FF_EAM) { - timer[FORCE] = computeForceEam(&eam, ¶m, &atom, &neighbor, &stats); - } else { - if(param.half_neigh) { - timer[FORCE] = computeForceLJHalfNeigh(¶m, &atom, &neighbor, &stats); - } else { - timer[FORCE] = computeForceLJFullNeigh(¶m, &atom, &neighbor, &stats); - } - } + timer[FORCE] = computeForce(&eam, ¶m, &atom, &neighbor, &stats); timer[NEIGH] = 0.0; timer[TOTAL] = getTimeStamp(); @@ -250,7 +259,6 @@ int main(int argc, char** argv) { for(int n = 0; n < param.ntimes; n++) { initialIntegrate(¶m, &atom); - if((n + 1) % param.reneigh_every) { updatePbc(&atom, ¶m); } else { @@ -261,16 +269,7 @@ int main(int argc, char** argv) { traceAddresses(¶m, &atom, &neighbor, n + 1); #endif - if(param.force_field == FF_EAM) { - timer[FORCE] += computeForceEam(&eam, ¶m, &atom, &neighbor, &stats); - } else { - if(param.half_neigh) { - timer[FORCE] += computeForceLJHalfNeigh(¶m, &atom, &neighbor, &stats); - } else { - timer[FORCE] += computeForceLJFullNeigh(¶m, &atom, &neighbor, &stats); - } - } - + timer[FORCE] += computeForce(&eam, ¶m, &atom, &neighbor, &stats); finalIntegrate(¶m, &atom); if(!((n + 1) % param.nstat) && (n+1) < param.ntimes) { diff --git a/lammps/util.c b/lammps/util.c index 0347ab6..a41c5a9 100644 --- a/lammps/util.c +++ b/lammps/util.c @@ -84,6 +84,7 @@ int str2ff(const char *string) { if(strncmp(string, "lj", 2) == 0) return FF_LJ; if(strncmp(string, "eam", 3) == 0) return FF_EAM; + if(strncmp(string, "dem", 3) == 0) return FF_DEM; return -1; } @@ -91,5 +92,6 @@ const char* ff2str(int ff) { if(ff == FF_LJ) { return "lj"; } if(ff == FF_EAM) { return "eam"; } + if(ff == FF_DEM) { return "dem"; } return "invalid"; }