moebiusband/MD-Bench
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Fix a few more bugs on gromacs variant

Browse Source 
Signed-off-by: Rafael Ravedutti <rafaelravedutti@gmail.com>
This commit is contained in:
Rafael Ravedutti 2022-01-31 23:46:20 +01:00
parent e0e6b6a68c
commit e64c3345bc
7 changed files with 85 additions and 28 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

11
gromacs/force_lj.c
View File

@ -67,7 +67,7 @@ double computeForceLJ(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *s
MD_FLOAT fiy = 0; MD_FLOAT fiy = 0;
MD_FLOAT fiz = 0; MD_FLOAT fiz = 0;
for(int cjj = 0; cjj < CLUSTER_DIM_N; cjj++) { for(int cjj = 0; cjj < atom->clusters[cj].natoms; cjj++) {
if(ci != cj || cii != cjj) { if(ci != cj || cii != cjj) {
MD_FLOAT delx = xtmp - cluster_x(cjptr, cjj); MD_FLOAT delx = xtmp - cluster_x(cjptr, cjj);
MD_FLOAT dely = ytmp - cluster_y(cjptr, cjj); MD_FLOAT dely = ytmp - cluster_y(cjptr, cjj);
@ -80,6 +80,11 @@ double computeForceLJ(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *s
fix += delx * force; fix += delx * force;
fiy += dely * force; fiy += dely * force;
fiz += delz * force; fiz += delz * force;
/*
if(force < -50.0 || force > 50.0) {
fprintf(stdout, "%d-%d/%d-%d: %f, %f, %f ---- %f, %f, %f, %f\n", ci, cii, cj, cjj, xtmp, ytmp, ztmp, fix, fiy, fiz, force);
}
*/
} }
} }
} }
@ -87,6 +92,10 @@ double computeForceLJ(Parameter *param, Atom *atom, Neighbor *neighbor, Stats *s
cluster_x(cifptr, cii) += fix; cluster_x(cifptr, cii) += fix;
cluster_y(cifptr, cii) += fiy; cluster_y(cifptr, cii) += fiy;
cluster_z(cifptr, cii) += fiz; cluster_z(cifptr, cii) += fiz;
/*
if(fix < -100.0 || fix > 100.0 || fiy < -100.0 || fiy > 100.0 || fiz < -100.0 || fiz > 100.0) {
fprintf(stdout, "%d-%d: %f, %f, %f ---- %f, %f, %f\n", ci, cii, xtmp, ytmp, ztmp, fix, fiy, fiz);
}*/
} }
} }

2
gromacs/includes/atom.h
View File

@ -25,8 +25,8 @@
#ifndef __ATOM_H_ #ifndef __ATOM_H_
#define __ATOM_H_ #define __ATOM_H_
#define CLUSTER_DIM_N 4
#define CLUSTER_DIM_M 4 #define CLUSTER_DIM_M 4
#define CLUSTER_DIM_N 4
typedef struct { typedef struct {
int bin; int bin;

3
gromacs/includes/vtk.h
View File

@ -24,5 +24,6 @@
#ifndef __VTK_H_ #ifndef __VTK_H_
#define __VTK_H_ #define __VTK_H_
extern int write_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep); extern int write_local_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep);
extern int write_ghost_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep);
#endif #endif

6
gromacs/main.c
View File

@ -260,7 +260,8 @@ int main(int argc, char** argv) {
timer[TOTAL] = getTimeStamp(); timer[TOTAL] = getTimeStamp();
if(param.vtk_file != NULL) { if(param.vtk_file != NULL) {
write_atoms_to_vtk_file(param.vtk_file, &atom, 0); write_local_atoms_to_vtk_file(param.vtk_file, &atom, 0);
write_ghost_atoms_to_vtk_file(param.vtk_file, &atom, 0);
} }
for(int n = 0; n < param.ntimes; n++) { for(int n = 0; n < param.ntimes; n++) {
@ -289,7 +290,8 @@ int main(int argc, char** argv) {
} }
if(param.vtk_file != NULL) { if(param.vtk_file != NULL) {
write_atoms_to_vtk_file(param.vtk_file, &atom, n + 1); write_local_atoms_to_vtk_file(param.vtk_file, &atom, n + 1);
write_ghost_atoms_to_vtk_file(param.vtk_file, &atom, n + 1);
} }
} }

19
gromacs/neighbor.c
View File

@ -289,7 +289,7 @@ void buildNeighbor(Atom *atom, Neighbor *neighbor) {
/* /*
printf("\ncutneighsq = %f, rbb_sq = %f\n", cutneighsq, rbb_sq); printf("\ncutneighsq = %f, rbb_sq = %f\n", cutneighsq, rbb_sq);
for(int ci = 4; ci < 6; ci++) { for(int ci = 0; ci < 6; ci++) {
MD_FLOAT *ciptr = cluster_pos_ptr(ci); MD_FLOAT *ciptr = cluster_pos_ptr(ci);
int* neighptr = &(neighbor->neighbors[ci * neighbor->maxneighs]); int* neighptr = &(neighbor->neighbors[ci * neighbor->maxneighs]);
printf("Cluster %d, bbx = {%f, %f}, bby = {%f, %f}, bbz = {%f, %f}\n", ci, atom->clusters[ci].bbminx, atom->clusters[ci].bbmaxx, atom->clusters[ci].bbminy, atom->clusters[ci].bbmaxy, atom->clusters[ci].bbminz, atom->clusters[ci].bbmaxz); printf("Cluster %d, bbx = {%f, %f}, bby = {%f, %f}, bbz = {%f, %f}\n", ci, atom->clusters[ci].bbminx, atom->clusters[ci].bbmaxx, atom->clusters[ci].bbminy, atom->clusters[ci].bbmaxy, atom->clusters[ci].bbminz, atom->clusters[ci].bbmaxz);
@ -377,7 +377,6 @@ void coord2bin2D(MD_FLOAT xin, MD_FLOAT yin, int *ix, int *iy) {
void binAtoms(Atom *atom) { void binAtoms(Atom *atom) {
printf("binAtoms start\n"); printf("binAtoms start\n");
int nall = atom->Nlocal + atom->Nghost;
int resize = 1; int resize = 1;
while(resize > 0) { while(resize > 0) {
@ -387,7 +386,7 @@ void binAtoms(Atom *atom) {
bincount[i] = 0; bincount[i] = 0;
} }
for(int i = 0; i < nall; i++) { for(int i = 0; i < atom->Nlocal; i++) {
int ibin = coord2bin(atom_x(i), atom_y(i)); int ibin = coord2bin(atom_x(i), atom_y(i));
if(ibin < 0 || ibin > mbins) { fprintf(stderr, "%d: %f, %f\n", i, atom_x(i), atom_y(i)); } if(ibin < 0 || ibin > mbins) { fprintf(stderr, "%d: %f, %f\n", i, atom_x(i), atom_y(i)); }
if(bincount[ibin] < atoms_per_bin) { if(bincount[ibin] < atoms_per_bin) {
@ -507,8 +506,15 @@ void buildClusters(Atom *atom) {
} }
} }
printf("buildClusters end\n");
}
void binClusters(Atom *atom) {
printf("binClusters start\n");
/* /*
for(int ci = 4; ci < 9; ci++) { fprintf(stdout, "Nghost = %d\n", atom->Nclusters_ghost);
for(int ci = atom->Nclusters_local; ci < atom->Nclusters_local + 4; ci++) {
MD_FLOAT *cptr = cluster_pos_ptr(ci); MD_FLOAT *cptr = cluster_pos_ptr(ci);
fprintf(stdout, "Cluster %d:\n", ci); fprintf(stdout, "Cluster %d:\n", ci);
fprintf(stdout, "bin=%d, Natoms=%d, bbox={%f,%f},{%f,%f},{%f,%f}\n", atom->clusters[ci].bin, atom->clusters[ci].natoms, atom->clusters[ci].bbminx, atom->clusters[ci].bbmaxx, atom->clusters[ci].bbminy, atom->clusters[ci].bbmaxy, atom->clusters[ci].bbminz, atom->clusters[ci].bbmaxz); fprintf(stdout, "bin=%d, Natoms=%d, bbox={%f,%f},{%f,%f},{%f,%f}\n", atom->clusters[ci].bin, atom->clusters[ci].natoms, atom->clusters[ci].bbminx, atom->clusters[ci].bbmaxx, atom->clusters[ci].bbminy, atom->clusters[ci].bbmaxy, atom->clusters[ci].bbminz, atom->clusters[ci].bbmaxz);
@ -518,11 +524,6 @@ void buildClusters(Atom *atom) {
} }
*/ */
printf("buildClusters end\n");
}
void binClusters(Atom *atom) {
printf("binClusters start\n");
const int nlocal = atom->Nclusters_local; const int nlocal = atom->Nclusters_local;
int resize = 1; int resize = 1;

25
gromacs/pbc.c
View File

@ -122,19 +122,18 @@ void updateAtomsPbc(Atom *atom, Parameter *param) {
* defining ghost atoms around domain * defining ghost atoms around domain
* only creates mapping and coordinate corrections * only creates mapping and coordinate corrections
* that are then enforced in updatePbc */ * that are then enforced in updatePbc */
#define ADDGHOST(dx,dy,dz) \ #define ADDGHOST(dx,dy,dz); \
Nghost++; \ Nghost++; \
border_map[Nghost] = ci; \ const int g_atom_idx = atom->Nclusters_local + Nghost; \
atom->PBCx[Nghost] = dx; \ border_map[Nghost] = ci; \
atom->PBCy[Nghost] = dy; \ atom->PBCx[Nghost] = dx; \
atom->PBCz[Nghost] = dz; \ atom->PBCy[Nghost] = dy; \
copy_cluster_types(atom, atom->Nclusters_local + Nghost, ci) atom->PBCz[Nghost] = dz; \
atom->clusters[g_atom_idx].natoms = atom->clusters[ci].natoms; \
void copy_cluster_types(Atom *atom, int dest, int src) { Nghost_atoms += atom->clusters[g_atom_idx].natoms; \
for(int cii = 0; cii < atom->clusters[src].natoms; cii++) { for(int cii = 0; cii < atom->clusters[ci].natoms; cii++) { \
atom->clusters[dest].type[cii] = atom->clusters[src].type[cii]; atom->clusters[g_atom_idx].type[cii] = atom->clusters[ci].type[cii]; \
} }
}
/* internal subroutines */ /* internal subroutines */
void growPbc(Atom* atom) { void growPbc(Atom* atom) {
@ -154,6 +153,7 @@ void setupPbc(Atom *atom, Parameter *param) {
MD_FLOAT zprd = param->zprd; MD_FLOAT zprd = param->zprd;
MD_FLOAT Cutneigh = param->cutneigh; MD_FLOAT Cutneigh = param->cutneigh;
int Nghost = -1; int Nghost = -1;
int Nghost_atoms = 0;
for(int ci = 0; ci < atom->Nclusters_local; ci++) { for(int ci = 0; ci < atom->Nclusters_local; ci++) {
if (atom->Nclusters_local + Nghost + 7 >= atom->Nclusters_max) { if (atom->Nclusters_local + Nghost + 7 >= atom->Nclusters_max) {
@ -205,6 +205,7 @@ void setupPbc(Atom *atom, Parameter *param) {
} }
// increase by one to make it the ghost atom count // increase by one to make it the ghost atom count
atom->Nghost = Nghost_atoms;
atom->Nclusters_ghost = Nghost + 1; atom->Nclusters_ghost = Nghost + 1;
atom->Nclusters = atom->Nclusters_local + Nghost + 1; atom->Nclusters = atom->Nclusters_local + Nghost + 1;

47
gromacs/vtk.c
View File

@ -3,9 +3,9 @@
#include <atom.h> #include <atom.h>
int write_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep) { int write_local_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep) {
char timestep_filename[128]; char timestep_filename[128];
snprintf(timestep_filename, sizeof timestep_filename, "%s_%d.vtk", filename, timestep); snprintf(timestep_filename, sizeof timestep_filename, "%s_local_%d.vtk", filename, timestep);
FILE* fp = fopen(timestep_filename, "wb"); FILE* fp = fopen(timestep_filename, "wb");
if(fp == NULL) { if(fp == NULL) {
@ -45,3 +45,46 @@ int write_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep) {
fclose(fp); fclose(fp);
return 0; return 0;
} }
int write_ghost_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep) {
char timestep_filename[128];
snprintf(timestep_filename, sizeof timestep_filename, "%s_ghost_%d.vtk", filename, timestep);
FILE* fp = fopen(timestep_filename, "wb");
if(fp == NULL) {
fprintf(stderr, "Could not open VTK file for writing!\n");
return -1;
}
fprintf(fp, "# vtk DataFile Version 2.0\n");
fprintf(fp, "Particle data\n");
fprintf(fp, "ASCII\n");
fprintf(fp, "DATASET UNSTRUCTURED_GRID\n");
fprintf(fp, "POINTS %d double\n", atom->Nghost);
for(int ci = atom->Nclusters_local; ci < atom->Nclusters_local + atom->Nclusters_ghost; ++ci) {
MD_FLOAT *cptr = cluster_pos_ptr(ci);
for(int cii = 0; cii < atom->clusters[ci].natoms; ++cii) {
fprintf(fp, "%.4f %.4f %.4f\n", cluster_x(cptr, cii), cluster_y(cptr, cii), cluster_z(cptr, cii));
}
}
fprintf(fp, "\n\n");
fprintf(fp, "CELLS %d %d\n", atom->Nghost, atom->Nghost * 2);
for(int i = 0; i < atom->Nghost; ++i) {
fprintf(fp, "1 %d\n", i);
}
fprintf(fp, "\n\n");
fprintf(fp, "CELL_TYPES %d\n", atom->Nghost);
for(int i = 0; i < atom->Nghost; ++i) {
fprintf(fp, "1\n");
}
fprintf(fp, "\n\n");
fprintf(fp, "POINT_DATA %d\n", atom->Nghost);
fprintf(fp, "SCALARS mass double\n");
fprintf(fp, "LOOKUP_TABLE default\n");
for(int i = 0; i < atom->Nghost; i++) {
fprintf(fp, "1.0\n");
}
fprintf(fp, "\n\n");
fclose(fp);
return 0;
}