Move common modules to common directory

Signed-off-by: Rafael Ravedutti <rafaelravedutti@gmail.com>
This commit is contained in:
Rafael Ravedutti 2022-08-17 17:56:31 +02:00
parent 47db9e86b0
commit 28d3946072
31 changed files with 40 additions and 1102 deletions

View File

@ -34,6 +34,7 @@ typedef struct {
char* param_file;
char* input_file;
char* vtk_file;
char* xtc_file;
MD_FLOAT epsilon;
MD_FLOAT sigma;
MD_FLOAT sigma6;
@ -44,6 +45,7 @@ typedef struct {
int ntimes;
int nstat;
int reneigh_every;
int prune_every;
int x_out_every;
int v_out_every;
int half_neigh;

View File

@ -31,6 +31,7 @@
void initParameter(Parameter *param) {
param->input_file = NULL;
param->vtk_file = NULL;
param->xtc_file = NULL;
param->eam_file = NULL;
param->force_field = FF_LJ;
param->epsilon = 1.0;
@ -54,6 +55,7 @@ void initParameter(Parameter *param) {
param->mass = 1.0;
param->dtforce = 0.5 * param->dt;
param->reneigh_every = 20;
param->prune_every = 1000;
param->x_out_every = 20;
param->v_out_every = 5;
param->half_neigh = 0;
@ -98,6 +100,7 @@ void readParameter(Parameter *param, const char *filename) {
PARSE_STRING(input_file);
PARSE_STRING(eam_file);
PARSE_STRING(vtk_file);
PARSE_STRING(xtc_file);
PARSE_REAL(epsilon);
PARSE_REAL(sigma);
PARSE_REAL(k_s);
@ -125,6 +128,7 @@ void readParameter(Parameter *param, const char *filename) {
PARSE_INT(pbc_z);
PARSE_INT(nstat);
PARSE_INT(reneigh_every);
PARSE_INT(prune_every);
PARSE_INT(x_out_every);
PARSE_INT(v_out_every);
PARSE_INT(half_neigh);
@ -150,12 +154,20 @@ void printParameter(Parameter *param) {
printf("VTK file: %s\n", param->vtk_file);
}
if(param->xtc_file != NULL) {
printf("XTC file: %s\n", param->xtc_file);
}
if(param->eam_file != NULL) {
printf("EAM file: %s\n", param->eam_file);
}
printf("\tForce field: %s\n", ff2str(param->force_field));
#ifdef CLUSTER_M
printf("\tKernel: %s, MxN: %dx%d, Vector width: %d\n", KERNEL_NAME, CLUSTER_M, CLUSTER_N, VECTOR_WIDTH);
#else
printf("\tKernel: %s\n", KERNEL_NAME);
#endif
printf("\tData layout: %s\n", POS_DATA_LAYOUT);
printf("\tFloating-point precision: %s\n", PRECISION_STRING);
printf("\tUnit cells (nx, ny, nz): %d, %d, %d\n", param->nx, param->ny, param->nz);
@ -173,6 +185,7 @@ void printParameter(Parameter *param) {
printf("\tNumber of timesteps: %d\n", param->ntimes);
printf("\tReport stats every (timesteps): %d\n", param->nstat);
printf("\tReneighbor every (timesteps): %d\n", param->reneigh_every);
printf("\tPrune every (timesteps): %d\n", param->prune_every);
printf("\tOutput positions every (timesteps): %d\n", param->x_out_every);
printf("\tOutput velocities every (timesteps): %d\n", param->v_out_every);
printf("\tDelta time (dt): %e\n", param->dt);

View File

@ -1,70 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <errno.h>
void* allocate (int alignment, size_t bytesize)
{
int errorCode;
void* ptr;
errorCode = posix_memalign(&ptr, alignment, bytesize);
if (errorCode) {
if (errorCode == EINVAL) {
fprintf(stderr,
"Error: Alignment parameter is not a power of two\n");
exit(EXIT_FAILURE);
}
if (errorCode == ENOMEM) {
fprintf(stderr,
"Error: Insufficient memory to fulfill the request\n");
exit(EXIT_FAILURE);
}
}
if (ptr == NULL) {
fprintf(stderr, "Error: posix_memalign failed!\n");
exit(EXIT_FAILURE);
}
return ptr;
}
void* reallocate (
void* ptr,
int alignment,
size_t newBytesize,
size_t oldBytesize)
{
void* newarray = allocate(alignment, newBytesize);
if(ptr != NULL) {
memcpy(newarray, ptr, oldBytesize);
free(ptr);
}
return newarray;
}

View File

@ -130,15 +130,31 @@ extern void growAtom(Atom*);
extern void growClusters(Atom*);
#ifdef AOS
#define POS_DATA_LAYOUT "AoS"
#define atom_x(i) atom->x[(i) * 3 + 0]
#define atom_y(i) atom->x[(i) * 3 + 1]
#define atom_z(i) atom->x[(i) * 3 + 2]
# define POS_DATA_LAYOUT "AoS"
# define atom_x(i) atom->x[(i) * 3 + 0]
# define atom_y(i) atom->x[(i) * 3 + 1]
# define atom_z(i) atom->x[(i) * 3 + 2]
/*
# define atom_vx(i) atom->vx[(i) * 3 + 0]
# define atom_vy(i) atom->vx[(i) * 3 + 1]
# define atom_vz(i) atom->vx[(i) * 3 + 2]
# define atom_fx(i) atom->fx[(i) * 3 + 0]
# define atom_fy(i) atom->fx[(i) * 3 + 1]
# define atom_fz(i) atom->fx[(i) * 3 + 2]
*/
#else
#define POS_DATA_LAYOUT "SoA"
#define atom_x(i) atom->x[i]
#define atom_y(i) atom->y[i]
#define atom_z(i) atom->z[i]
# define POS_DATA_LAYOUT "SoA"
# define atom_x(i) atom->x[i]
# define atom_y(i) atom->y[i]
# define atom_z(i) atom->z[i]
#endif
// TODO: allow to switch velocites and forces to AoS
# define atom_vx(i) atom->vx[i]
# define atom_vy(i) atom->vy[i]
# define atom_vz(i) atom->vz[i]
# define atom_fx(i) atom->fx[i]
# define atom_fy(i) atom->fy[i]
# define atom_fz(i) atom->fz[i]
#endif

View File

@ -1,69 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#ifndef __PARAMETER_H_
#define __PARAMETER_H_
#if PRECISION == 1
#define MD_FLOAT float
#else
#define MD_FLOAT double
#endif
typedef struct {
int force_field;
char* param_file;
char* input_file;
char* vtk_file;
char* xtc_file;
MD_FLOAT epsilon;
MD_FLOAT sigma;
MD_FLOAT sigma6;
MD_FLOAT temp;
MD_FLOAT rho;
MD_FLOAT mass;
int ntypes;
int ntimes;
int nstat;
int reneigh_every;
int prune_every;
int x_out_every;
int v_out_every;
int half_neigh;
int nx, ny, nz;
MD_FLOAT dt;
MD_FLOAT dtforce;
MD_FLOAT cutforce;
MD_FLOAT skin;
MD_FLOAT cutneigh;
MD_FLOAT lattice;
MD_FLOAT xlo, xhi, ylo, yhi, zlo, zhi;
MD_FLOAT xprd, yprd, zprd;
double proc_freq;
char* eam_file;
} Parameter;
void initParameter(Parameter*);
void readParameter(Parameter*, const char*);
void printParameter(Parameter*);
#endif

View File

@ -1,163 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
//---
#include <atom.h>
#include <parameter.h>
#include <util.h>
void initParameter(Parameter *param) {
param->input_file = NULL;
param->vtk_file = NULL;
param->xtc_file = NULL;
param->eam_file = NULL;
param->force_field = FF_LJ;
param->epsilon = 1.0;
param->sigma = 1.0;
param->sigma6 = 1.0;
param->rho = 0.8442;
param->ntypes = 4;
param->ntimes = 200;
param->dt = 0.005;
param->nx = 32;
param->ny = 32;
param->nz = 32;
param->half_neigh = 0;
param->cutforce = 2.5;
param->skin = 0.3;
param->cutneigh = param->cutforce + param->skin;
param->temp = 1.44;
param->nstat = 100;
param->mass = 1.0;
param->dtforce = 0.5 * param->dt;
param->reneigh_every = 20;
param->prune_every = 1000;
param->x_out_every = 20;
param->v_out_every = 5;
param->proc_freq = 2.4;
}
void readParameter(Parameter *param, const char *filename) {
FILE *fp = fopen(filename, "r");
char line[MAXLINE];
int i;
if(!fp) {
fprintf(stderr, "Could not open parameter file: %s\n", filename);
exit(-1);
}
while(!feof(fp)) {
line[0] = '\0';
fgets(line, MAXLINE, fp);
for(i = 0; line[i] != '\0' && line[i] != '#'; i++);
line[i] = '\0';
char *tok = strtok(line, " ");
char *val = strtok(NULL, " ");
#define PARSE_PARAM(p,f) if(strncmp(tok, #p, sizeof(#p) / sizeof(#p[0]) - 1) == 0) { param->p = f(val); }
#define PARSE_STRING(p) PARSE_PARAM(p, strdup)
#define PARSE_INT(p) PARSE_PARAM(p, atoi)
#define PARSE_REAL(p) PARSE_PARAM(p, atof)
if(tok != NULL && val != NULL) {
PARSE_PARAM(force_field, str2ff);
PARSE_STRING(input_file);
PARSE_STRING(eam_file);
PARSE_STRING(vtk_file);
PARSE_STRING(xtc_file);
PARSE_REAL(epsilon);
PARSE_REAL(sigma);
PARSE_REAL(rho);
PARSE_REAL(dt);
PARSE_REAL(cutforce);
PARSE_REAL(skin);
PARSE_REAL(temp);
PARSE_REAL(mass);
PARSE_REAL(proc_freq);
PARSE_INT(ntypes);
PARSE_INT(ntimes);
PARSE_INT(nx);
PARSE_INT(ny);
PARSE_INT(nz);
PARSE_INT(half_neigh);
PARSE_INT(nstat);
PARSE_INT(reneigh_every);
PARSE_INT(prune_every);
PARSE_INT(x_out_every);
PARSE_INT(v_out_every);
}
}
// Update sigma6 parameter
MD_FLOAT s2 = param->sigma * param->sigma;
param->sigma6 = s2 * s2 * s2;
fclose(fp);
}
void printParameter(Parameter *param) {
printf("Parameters:\n");
if(param->input_file != NULL) {
printf("Input file: %s\n", param->input_file);
}
if(param->vtk_file != NULL) {
printf("VTK file: %s\n", param->vtk_file);
}
if(param->xtc_file != NULL) {
printf("XTC file: %s\n", param->xtc_file);
}
if(param->eam_file != NULL) {
printf("EAM file: %s\n", param->eam_file);
}
printf("\tForce field: %s\n", ff2str(param->force_field));
printf("\tKernel: %s, MxN: %dx%d, Vector width: %d\n", KERNEL_NAME, CLUSTER_M, CLUSTER_N, VECTOR_WIDTH);
printf("\tData layout: %s\n", POS_DATA_LAYOUT);
printf("\tFloating-point precision: %s\n", PRECISION_STRING);
printf("\tUnit cells (nx, ny, nz): %d, %d, %d\n", param->nx, param->ny, param->nz);
printf("\tDomain box sizes (x, y, z): %e, %e, %e\n", param->xprd, param->yprd, param->zprd);
printf("\tLattice size: %e\n", param->lattice);
printf("\tEpsilon: %e\n", param->epsilon);
printf("\tSigma: %e\n", param->sigma);
printf("\tTemperature: %e\n", param->temp);
printf("\tRHO: %e\n", param->rho);
printf("\tMass: %e\n", param->mass);
printf("\tNumber of types: %d\n", param->ntypes);
printf("\tNumber of timesteps: %d\n", param->ntimes);
printf("\tReport stats every (timesteps): %d\n", param->nstat);
printf("\tReneighbor every (timesteps): %d\n", param->reneigh_every);
printf("\tPrune every (timesteps): %d\n", param->prune_every);
printf("\tOutput positions every (timesteps): %d\n", param->x_out_every);
printf("\tOutput velocities every (timesteps): %d\n", param->v_out_every);
printf("\tDelta time (dt): %e\n", param->dt);
printf("\tCutoff radius: %e\n", param->cutforce);
printf("\tSkin: %e\n", param->skin);
printf("\tHalf neighbor lists: %d\n", param->half_neigh);
printf("\tProcessor frequency (GHz): %.4f\n", param->proc_freq);
}

View File

@ -1,136 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <thermo.h>
#include <util.h>
static int *steparr;
static MD_FLOAT *tmparr;
static MD_FLOAT *engarr;
static MD_FLOAT *prsarr;
static MD_FLOAT mvv2e;
static MD_FLOAT dof_boltz;
static MD_FLOAT t_scale;
static MD_FLOAT p_scale;
static MD_FLOAT e_scale;
static MD_FLOAT t_act;
static MD_FLOAT p_act;
static MD_FLOAT e_act;
static int mstat;
/* exported subroutines */
void setupThermo(Parameter *param, int natoms)
{
int maxstat = param->ntimes / param->nstat + 2;
steparr = (int*) malloc(maxstat * sizeof(int));
tmparr = (MD_FLOAT*) malloc(maxstat * sizeof(MD_FLOAT));
engarr = (MD_FLOAT*) malloc(maxstat * sizeof(MD_FLOAT));
prsarr = (MD_FLOAT*) malloc(maxstat * sizeof(MD_FLOAT));
if(param->force_field == FF_LJ) {
mvv2e = 1.0;
dof_boltz = (natoms * 3 - 3);
t_scale = mvv2e / dof_boltz;
p_scale = 1.0 / 3 / param->xprd / param->yprd / param->zprd;
e_scale = 0.5;
} else {
mvv2e = 1.036427e-04;
dof_boltz = (natoms * 3 - 3) * 8.617343e-05;
t_scale = mvv2e / dof_boltz;
p_scale = 1.602176e+06 / 3 / param->xprd / param->yprd / param->zprd;
e_scale = 524287.985533;//16.0;
param->dtforce /= mvv2e;
}
}
void computeThermo(int iflag, Parameter *param, Atom *atom)
{
MD_FLOAT t = 0.0, p;
MD_FLOAT* vx = atom->vx;
MD_FLOAT* vy = atom->vy;
MD_FLOAT* vz = atom->vz;
for(int i = 0; i < atom->Nlocal; i++) {
t += (vx[i] * vx[i] + vy[i] * vy[i] + vz[i] * vz[i]) * param->mass;
}
t = t * t_scale;
p = (t * dof_boltz) * p_scale;
int istep = iflag;
if(iflag == -1){
istep = param->ntimes;
}
if(iflag == 0){
mstat = 0;
}
steparr[mstat] = istep;
tmparr[mstat] = t;
prsarr[mstat] = p;
mstat++;
fprintf(stdout, "%i\t%e\t%e\n", istep, t, p);
}
void adjustThermo(Parameter *param, Atom *atom)
{
/* zero center-of-mass motion */
MD_FLOAT vxtot = 0.0; MD_FLOAT vytot = 0.0; MD_FLOAT vztot = 0.0;
MD_FLOAT* vx = atom->vx; MD_FLOAT* vy = atom->vy; MD_FLOAT* vz = atom->vz;
for(int i = 0; i < atom->Nlocal; i++) {
vxtot += vx[i];
vytot += vy[i];
vztot += vz[i];
}
vxtot = vxtot / atom->Natoms;
vytot = vytot / atom->Natoms;
vztot = vztot / atom->Natoms;
for(int i = 0; i < atom->Nlocal; i++) {
vx[i] -= vxtot;
vy[i] -= vytot;
vz[i] -= vztot;
}
t_act = 0;
MD_FLOAT t = 0.0;
for(int i = 0; i < atom->Nlocal; i++) {
t += (vx[i] * vx[i] + vy[i] * vy[i] + vz[i] * vz[i]) * param->mass;
}
t *= t_scale;
MD_FLOAT factor = sqrt(param->temp / t);
for(int i = 0; i < atom->Nlocal; i++) {
vx[i] *= factor;
vy[i] *= factor;
vz[i] *= factor;
}
}

View File

@ -1,285 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <allocate.h>
#include <atom.h>
#include <eam.h>
#include <parameter.h>
#include <util.h>
#ifndef MAXLINE
#define MAXLINE 4096
#endif
void initEam(Eam* eam, Parameter* param) {
int ntypes = param->ntypes;
eam->nmax = 0;
eam->fp = NULL;
coeff(eam, param);
init_style(eam, param);
}
void coeff(Eam* eam, Parameter* param) {
read_eam_file(&eam->file, param->eam_file);
param->mass = eam->file.mass;
param->cutforce = eam->file.cut;
param->cutneigh = param->cutforce + 1.0;
param->temp = 600.0;
param->dt = 0.001;
param->rho = 0.07041125;
param->dtforce = 0.5 * param->dt / param->mass;
}
void init_style(Eam* eam, Parameter* param) {
// convert read-in file(s) to arrays and spline them
file2array(eam);
array2spline(eam, param);
}
void read_eam_file(Funcfl* file, const char* filename) {
FILE* fptr;
char line[MAXLINE];
fptr = fopen(filename, "r");
if(fptr == NULL) {
printf("Can't open EAM Potential file: %s\n", filename);
exit(0);
}
int tmp;
fgets(line, MAXLINE, fptr);
fgets(line, MAXLINE, fptr);
sscanf(line, "%d %lg", &tmp, &(file->mass));
fgets(line, MAXLINE, fptr);
sscanf(line, "%d %lg %d %lg %lg", &file->nrho, &file->drho, &file->nr, &file->dr, &file->cut);
//printf("Read: %lf %i %lf %i %lf %lf\n",file->mass,file->nrho,file->drho,file->nr,file->dr,file->cut);
file->frho = (MD_FLOAT *) allocate(ALIGNMENT, (file->nrho + 1) * sizeof(MD_FLOAT));
file->rhor = (MD_FLOAT *) allocate(ALIGNMENT, (file->nr + 1) * sizeof(MD_FLOAT));
file->zr = (MD_FLOAT *) allocate(ALIGNMENT, (file->nr + 1) * sizeof(MD_FLOAT));
grab(fptr, file->nrho, file->frho);
grab(fptr, file->nr, file->zr);
grab(fptr, file->nr, file->rhor);
for(int i = file->nrho; i > 0; i--) file->frho[i] = file->frho[i - 1];
for(int i = file->nr; i > 0; i--) file->rhor[i] = file->rhor[i - 1];
for(int i = file->nr; i > 0; i--) file->zr[i] = file->zr[i - 1];
fclose(fptr);
}
void file2array(Eam* eam) {
int i, j, k, m, n;
double sixth = 1.0 / 6.0;
// determine max function params from all active funcfl files
// active means some element is pointing at it via map
int active;
double rmax, rhomax;
eam->dr = eam->drho = rmax = rhomax = 0.0;
active = 0;
Funcfl* file = &eam->file;
eam->dr = MAX(eam->dr, file->dr);
eam->drho = MAX(eam->drho, file->drho);
rmax = MAX(rmax, (file->nr - 1) * file->dr);
rhomax = MAX(rhomax, (file->nrho - 1) * file->drho);
// set nr,nrho from cutoff and spacings
// 0.5 is for round-off in divide
eam->nr = (int)(rmax / eam->dr + 0.5);
eam->nrho = (int)(rhomax / eam->drho + 0.5);
// ------------------------------------------------------------------
// setup frho arrays
// ------------------------------------------------------------------
// allocate frho arrays
// nfrho = # of funcfl files + 1 for zero array
eam->frho = (MD_FLOAT *) allocate(ALIGNMENT, (eam->nrho + 1) * sizeof(MD_FLOAT));
// interpolate each file's frho to a single grid and cutoff
double r, p, cof1, cof2, cof3, cof4;
n = 0;
for(m = 1; m <= eam->nrho; m++) {
r = (m - 1) * eam->drho;
p = r / file->drho + 1.0;
k = (int)(p);
k = MIN(k, file->nrho - 2);
k = MAX(k, 2);
p -= k;
p = MIN(p, 2.0);
cof1 = -sixth * p * (p - 1.0) * (p - 2.0);
cof2 = 0.5 * (p * p - 1.0) * (p - 2.0);
cof3 = -0.5 * p * (p + 1.0) * (p - 2.0);
cof4 = sixth * p * (p * p - 1.0);
eam->frho[m] = cof1 * file->frho[k - 1] + cof2 * file->frho[k] +
cof3 * file->frho[k + 1] + cof4 * file->frho[k + 2];
}
// ------------------------------------------------------------------
// setup rhor arrays
// ------------------------------------------------------------------
// allocate rhor arrays
// nrhor = # of funcfl files
eam->rhor = (MD_FLOAT *) allocate(ALIGNMENT, (eam->nr + 1) * sizeof(MD_FLOAT));
// interpolate each file's rhor to a single grid and cutoff
for(m = 1; m <= eam->nr; m++) {
r = (m - 1) * eam->dr;
p = r / file->dr + 1.0;
k = (int)(p);
k = MIN(k, file->nr - 2);
k = MAX(k, 2);
p -= k;
p = MIN(p, 2.0);
cof1 = -sixth * p * (p - 1.0) * (p - 2.0);
cof2 = 0.5 * (p * p - 1.0) * (p - 2.0);
cof3 = -0.5 * p * (p + 1.0) * (p - 2.0);
cof4 = sixth * p * (p * p - 1.0);
eam->rhor[m] = cof1 * file->rhor[k - 1] + cof2 * file->rhor[k] +
cof3 * file->rhor[k + 1] + cof4 * file->rhor[k + 2];
//if(m==119)printf("BuildRho: %e %e %e %e %e %e\n",rhor[m],cof1,cof2,cof3,cof4,file->rhor[k]);
}
// type2rhor[i][j] = which rhor array (0 to nrhor-1) each type pair maps to
// for funcfl files, I,J mapping only depends on I
// OK if map = -1 (non-EAM atom in pair hybrid) b/c type2rhor not used
// ------------------------------------------------------------------
// setup z2r arrays
// ------------------------------------------------------------------
// allocate z2r arrays
// nz2r = N*(N+1)/2 where N = # of funcfl files
eam->z2r = (MD_FLOAT *) allocate(ALIGNMENT, (eam->nr + 1) * sizeof(MD_FLOAT));
// create a z2r array for each file against other files, only for I >= J
// interpolate zri and zrj to a single grid and cutoff
double zri, zrj;
Funcfl* ifile = &eam->file;
Funcfl* jfile = &eam->file;
for(m = 1; m <= eam->nr; m++) {
r = (m - 1) * eam->dr;
p = r / ifile->dr + 1.0;
k = (int)(p);
k = MIN(k, ifile->nr - 2);
k = MAX(k, 2);
p -= k;
p = MIN(p, 2.0);
cof1 = -sixth * p * (p - 1.0) * (p - 2.0);
cof2 = 0.5 * (p * p - 1.0) * (p - 2.0);
cof3 = -0.5 * p * (p + 1.0) * (p - 2.0);
cof4 = sixth * p * (p * p - 1.0);
zri = cof1 * ifile->zr[k - 1] + cof2 * ifile->zr[k] +
cof3 * ifile->zr[k + 1] + cof4 * ifile->zr[k + 2];
p = r / jfile->dr + 1.0;
k = (int)(p);
k = MIN(k, jfile->nr - 2);
k = MAX(k, 2);
p -= k;
p = MIN(p, 2.0);
cof1 = -sixth * p * (p - 1.0) * (p - 2.0);
cof2 = 0.5 * (p * p - 1.0) * (p - 2.0);
cof3 = -0.5 * p * (p + 1.0) * (p - 2.0);
cof4 = sixth * p * (p * p - 1.0);
zrj = cof1 * jfile->zr[k - 1] + cof2 * jfile->zr[k] +
cof3 * jfile->zr[k + 1] + cof4 * jfile->zr[k + 2];
eam->z2r[m] = 27.2 * 0.529 * zri * zrj;
}
}
void array2spline(Eam* eam, Parameter* param) {
eam->rdr = 1.0 / eam->dr;
eam->rdrho = 1.0 / eam->drho;
eam->nrho_tot = (eam->nrho + 1) * 7 + 64;
eam->nr_tot = (eam->nr + 1) * 7 + 64;
eam->nrho_tot -= eam->nrho_tot%64;
eam->nr_tot -= eam->nr_tot%64;
int ntypes = param->ntypes;
eam->frho_spline = (MD_FLOAT *) allocate(ALIGNMENT, ntypes * ntypes * eam->nrho_tot * sizeof(MD_FLOAT));
eam->rhor_spline = (MD_FLOAT *) allocate(ALIGNMENT, ntypes * ntypes * eam->nr_tot * sizeof(MD_FLOAT));
eam->z2r_spline = (MD_FLOAT *) allocate(ALIGNMENT, ntypes * ntypes * eam->nr_tot * sizeof(MD_FLOAT));
interpolate(eam->nrho, eam->drho, eam->frho, eam->frho_spline);
interpolate(eam->nr, eam->dr, eam->rhor, eam->rhor_spline);
interpolate(eam->nr, eam->dr, eam->z2r, eam->z2r_spline);
// replicate data for multiple types;
for(int tt = 0; tt < ntypes * ntypes; tt++) {
for(int k = 0; k < eam->nrho_tot; k++)
eam->frho_spline[tt*eam->nrho_tot + k] = eam->frho_spline[k];
for(int k = 0; k < eam->nr_tot; k++)
eam->rhor_spline[tt*eam->nr_tot + k] = eam->rhor_spline[k];
for(int k = 0; k < eam->nr_tot; k++)
eam->z2r_spline[tt*eam->nr_tot + k] = eam->z2r_spline[k];
}
}
void interpolate(int n, MD_FLOAT delta, MD_FLOAT* f, MD_FLOAT* spline) {
for(int m = 1; m <= n; m++) spline[m * 7 + 6] = f[m];
spline[1 * 7 + 5] = spline[2 * 7 + 6] - spline[1 * 7 + 6];
spline[2 * 7 + 5] = 0.5 * (spline[3 * 7 + 6] - spline[1 * 7 + 6]);
spline[(n - 1) * 7 + 5] = 0.5 * (spline[n * 7 + 6] - spline[(n - 2) * 7 + 6]);
spline[n * 7 + 5] = spline[n * 7 + 6] - spline[(n - 1) * 7 + 6];
for(int m = 3; m <= n - 2; m++)
spline[m * 7 + 5] = ((spline[(m - 2) * 7 + 6] - spline[(m + 2) * 7 + 6]) +
8.0 * (spline[(m + 1) * 7 + 6] - spline[(m - 1) * 7 + 6])) / 12.0;
for(int m = 1; m <= n - 1; m++) {
spline[m * 7 + 4] = 3.0 * (spline[(m + 1) * 7 + 6] - spline[m * 7 + 6]) -
2.0 * spline[m * 7 + 5] - spline[(m + 1) * 7 + 5];
spline[m * 7 + 3] = spline[m * 7 + 5] + spline[(m + 1) * 7 + 5] -
2.0 * (spline[(m + 1) * 7 + 6] - spline[m * 7 + 6]);
}
spline[n * 7 + 4] = 0.0;
spline[n * 7 + 3] = 0.0;
for(int m = 1; m <= n; m++) {
spline[m * 7 + 2] = spline[m * 7 + 5] / delta;
spline[m * 7 + 1] = 2.0 * spline[m * 7 + 4] / delta;
spline[m * 7 + 0] = 3.0 * spline[m * 7 + 3] / delta;
}
}
void grab(FILE* fptr, int n, MD_FLOAT* list) {
char* ptr;
char line[MAXLINE];
int i = 0;
while(i < n) {
fgets(line, MAXLINE, fptr);
ptr = strtok(line, " \t\n\r\f");
list[i++] = atof(ptr);
while(ptr = strtok(NULL, " \t\n\r\f")) list[i++] = atof(ptr);
}
}

View File

@ -30,9 +30,8 @@
#include <stats.h>
#include <timing.h>
// TODO: Joint common files for gromacs and lammps variants
#ifdef SIMD_KERNEL_AVAILABLE
#include "../gromacs/includes/simd.h"
#include <simd.h>
#endif
double computeForceLJFullNeigh_plain_c(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *stats) {

View File

@ -1,29 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdlib.h>
#ifndef __ALLOCATE_H_
#define __ALLOCATE_H_
extern void* allocate (int alignment, size_t bytesize);
extern void* reallocate (void* ptr, int alignment, size_t newBytesize, size_t oldBytesize);
#endif

View File

@ -1,55 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdio.h>
#include <atom.h>
#include <parameter.h>
#ifndef __EAM_H_
#define __EAM_H_
typedef struct {
int nrho, nr;
MD_FLOAT drho, dr, cut, mass;
MD_FLOAT *frho, *rhor, *zr;
} Funcfl;
typedef struct {
MD_FLOAT* fp;
int nmax;
int nrho, nr;
int nrho_tot, nr_tot;
MD_FLOAT dr, rdr, drho, rdrho;
MD_FLOAT *frho, *rhor, *z2r;
MD_FLOAT *rhor_spline, *frho_spline, *z2r_spline;
Funcfl file;
} Eam;
void initEam(Eam* eam, Parameter* param);
void coeff(Eam* eam, Parameter* param);
void init_style(Eam* eam, Parameter *param);
void read_eam_file(Funcfl* file, const char* filename);
void file2array(Eam* eam);
void array2spline(Eam* eam, Parameter* param);
void interpolate(int n, MD_FLOAT delta, MD_FLOAT* f, MD_FLOAT* spline);
void grab(FILE* fptr, int n, MD_FLOAT* list);
#endif

View File

@ -1,170 +0,0 @@
/*
* =======================================================================================
*
* Filename: likwid-marker.h
*
* Description: Header File of likwid Marker API
*
* Version: <VERSION>
* Released: <DATE>
*
* Authors: Thomas Gruber (tg), thomas.roehl@googlemail.com
*
* Project: likwid
*
* Copyright (C) 2016 RRZE, University Erlangen-Nuremberg
*
* This program is free software: you can redistribute it and/or modify it under
* the terms of the GNU General Public License as published by the Free Software
* Foundation, either version 3 of the License, or (at your option) any later
* version.
*
* This program is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License along with
* this program. If not, see <http://www.gnu.org/licenses/>.
*
* =======================================================================================
*/
#ifndef LIKWID_MARKER_H
#define LIKWID_MARKER_H
/** \addtogroup MarkerAPI Marker API module
* @{
*/
/*!
\def LIKWID_MARKER_INIT
Shortcut for likwid_markerInit() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_THREADINIT
Shortcut for likwid_markerThreadInit() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_REGISTER(regionTag)
Shortcut for likwid_markerRegisterRegion() with \a regionTag if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_START(regionTag)
Shortcut for likwid_markerStartRegion() with \a regionTag if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_STOP(regionTag)
Shortcut for likwid_markerStopRegion() with \a regionTag if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_GET(regionTag, nevents, events, time, count)
Shortcut for likwid_markerGetResults() for \a regionTag if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_SWITCH
Shortcut for likwid_markerNextGroup() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_RESET(regionTag)
Shortcut for likwid_markerResetRegion() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_MARKER_CLOSE
Shortcut for likwid_markerClose() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/** @}*/
#ifdef LIKWID_PERFMON
#include <likwid.h>
#define LIKWID_MARKER_INIT likwid_markerInit()
#define LIKWID_MARKER_THREADINIT likwid_markerThreadInit()
#define LIKWID_MARKER_SWITCH likwid_markerNextGroup()
#define LIKWID_MARKER_REGISTER(regionTag) likwid_markerRegisterRegion(regionTag)
#define LIKWID_MARKER_START(regionTag) likwid_markerStartRegion(regionTag)
#define LIKWID_MARKER_STOP(regionTag) likwid_markerStopRegion(regionTag)
#define LIKWID_MARKER_CLOSE likwid_markerClose()
#define LIKWID_MARKER_RESET(regionTag) likwid_markerResetRegion(regionTag)
#define LIKWID_MARKER_GET(regionTag, nevents, events, time, count) likwid_markerGetRegion(regionTag, nevents, events, time, count)
#else /* LIKWID_PERFMON */
#define LIKWID_MARKER_INIT
#define LIKWID_MARKER_THREADINIT
#define LIKWID_MARKER_SWITCH
#define LIKWID_MARKER_REGISTER(regionTag)
#define LIKWID_MARKER_START(regionTag)
#define LIKWID_MARKER_STOP(regionTag)
#define LIKWID_MARKER_CLOSE
#define LIKWID_MARKER_GET(regionTag, nevents, events, time, count)
#define LIKWID_MARKER_RESET(regionTag)
#endif /* LIKWID_PERFMON */
/** \addtogroup NvMarkerAPI NvMarker API module (MarkerAPI for Nvidia GPUs)
* @{
*/
/*!
\def LIKWID_NVMARKER_INIT
Shortcut for likwid_gpuMarkerInit() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_THREADINIT
Shortcut for likwid_gpuMarkerThreadInit() if compiled with -DLIKWID_PERFMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_REGISTER(regionTag)
Shortcut for likwid_gpuMarkerRegisterRegion() with \a regionTag if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_START(regionTag)
Shortcut for likwid_gpuMarkerStartRegion() with \a regionTag if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_STOP(regionTag)
Shortcut for likwid_gpuMarkerStopRegion() with \a regionTag if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_GET(regionTag, ngpus, nevents, events, time, count)
Shortcut for likwid_gpuMarkerGetRegion() for \a regionTag if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_SWITCH
Shortcut for likwid_gpuMarkerNextGroup() if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_RESET(regionTag)
Shortcut for likwid_gpuMarkerResetRegion() if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/*!
\def LIKWID_NVMARKER_CLOSE
Shortcut for likwid_gpuMarkerClose() if compiled with -DLIKWID_NVMON. Otherwise no operation is performed
*/
/** @}*/
#ifdef LIKWID_NVMON
#ifndef LIKWID_WITH_NVMON
#define LIKWID_WITH_NVMON
#endif
#include <likwid.h>
#define LIKWID_NVMARKER_INIT likwid_gpuMarkerInit()
#define LIKWID_NVMARKER_THREADINIT likwid_gpuMarkerThreadInit()
#define LIKWID_NVMARKER_SWITCH likwid_gpuMarkerNextGroup()
#define LIKWID_NVMARKER_REGISTER(regionTag) likwid_gpuMarkerRegisterRegion(regionTag)
#define LIKWID_NVMARKER_START(regionTag) likwid_gpuMarkerStartRegion(regionTag)
#define LIKWID_NVMARKER_STOP(regionTag) likwid_gpuMarkerStopRegion(regionTag)
#define LIKWID_NVMARKER_CLOSE likwid_gpuMarkerClose()
#define LIKWID_NVMARKER_RESET(regionTag) likwid_gpuMarkerResetRegion(regionTag)
#define LIKWID_NVMARKER_GET(regionTag, ngpus, nevents, events, time, count) \
likwid_gpuMarkerGetRegion(regionTag, ngpus, nevents, events, time, count)
#else /* LIKWID_NVMON */
#define LIKWID_NVMARKER_INIT
#define LIKWID_NVMARKER_THREADINIT
#define LIKWID_NVMARKER_SWITCH
#define LIKWID_NVMARKER_REGISTER(regionTag)
#define LIKWID_NVMARKER_START(regionTag)
#define LIKWID_NVMARKER_STOP(regionTag)
#define LIKWID_NVMARKER_CLOSE
#define LIKWID_NVMARKER_GET(regionTag, nevents, events, time, count)
#define LIKWID_NVMARKER_RESET(regionTag)
#endif /* LIKWID_NVMON */
#endif /* LIKWID_MARKER_H */

View File

@ -1,31 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <parameter.h>
#include <atom.h>
#ifndef __THERMO_H_
#define __THERMO_H_
extern void setupThermo(Parameter*, int);
extern void computeThermo(int, Parameter*, Atom*);
extern void adjustThermo(Parameter*, Atom*);
#endif

View File

@ -1,11 +0,0 @@
#ifndef __TIMERS_H_
#define __TIMERS_H_
typedef enum {
TOTAL = 0,
NEIGH,
FORCE,
NUMTIMER
} timertype;
#endif

View File

@ -1,30 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#ifndef __TIMING_H_
#define __TIMING_H_
extern double getTimeStamp();
extern double getTimeResolution();
extern double getTimeStamp_();
#endif

View File

@ -1,43 +0,0 @@
/*
* =======================================================================================
*
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
*
* This file is part of MD-Bench.
*
* MD-Bench is free software: you can redistribute it and/or modify it
* under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* MD-Bench is distributed in the hope that it will be useful, but WITHOUT ANY
* WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
* PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License along
* with MD-Bench. If not, see <https://www.gnu.org/licenses/>.
* =======================================================================================
*/
#include <stdlib.h>
#include <time.h>
double getTimeStamp()
{
struct timespec ts;
clock_gettime(CLOCK_MONOTONIC, &ts);
return (double)ts.tv_sec + (double)ts.tv_nsec * 1.e-9;
}
double getTimeResolution()
{
struct timespec ts;
clock_getres(CLOCK_MONOTONIC, &ts);
return (double)ts.tv_sec + (double)ts.tv_nsec * 1.e-9;
}
double getTimeStamp_()
{
return getTimeStamp();
}