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

Merge pull request #2 from RRZE-HPC/aos

Browse Source 
Aos
This commit is contained in:
moebiusband73 2021-03-24 08:50:29 +01:00 committed by GitHub
commit dcfa1c4135
12 changed files with 124 additions and 111 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

15
Makefile
View File

@ -1,5 +1,3 @@
TAG = CLANG
#CONFIGURE BUILD SYSTEM #CONFIGURE BUILD SYSTEM
TARGET = MDBench-$(TAG) TARGET = MDBench-$(TAG)
BUILD_DIR = ./$(TAG) BUILD_DIR = ./$(TAG)
@ -8,14 +6,23 @@ MAKE_DIR = ./
Q ?= @ Q ?= @
#DO NOT EDIT BELOW #DO NOT EDIT BELOW
include $(MAKE_DIR)/config.mk
include $(MAKE_DIR)/include_$(TAG).mk include $(MAKE_DIR)/include_$(TAG).mk
INCLUDES += -I./src/includes INCLUDES += -I./src/includes
ifeq ($(strip $(DATA_LAYOUT)),AOS)
DEFINES += -DAOS
endif
ifeq ($(strip $(DATA_TYPE)),SP)
DEFINES += -DPRECISION=1
else
DEFINES += -DPRECISION=2
endif
VPATH = $(SRC_DIR) VPATH = $(SRC_DIR)
ASM = $(patsubst $(SRC_DIR)/%.c, $(BUILD_DIR)/%.s,$(wildcard $(SRC_DIR)/*.c)) ASM = $(patsubst $(SRC_DIR)/%.c, $(BUILD_DIR)/%.s,$(wildcard $(SRC_DIR)/*.c))
OBJ = $(patsubst $(SRC_DIR)/%.c, $(BUILD_DIR)/%.o,$(wildcard $(SRC_DIR)/*.c)) OBJ = $(patsubst $(SRC_DIR)/%.c, $(BUILD_DIR)/%.o,$(wildcard $(SRC_DIR)/*.c))
CPPFLAGS := $(CPPFLAGS) $(DEFINES) $(INCLUDES) CPPFLAGS := $(CPPFLAGS) $(DEFINES) $(OPTIONS) $(INCLUDES)
${TARGET}: $(BUILD_DIR) $(OBJ) ${TARGET}: $(BUILD_DIR) $(OBJ)
@echo "===> LINKING $(TARGET)" @echo "===> LINKING $(TARGET)"

7
README.md
View File

@ -4,9 +4,10 @@ A simple, sequential C implementation of the [Mantevo miniMD](https://github.co
## Build ## Build
1. Open the `Makefile` and edit the `TAG` value according to the tool chain used. Currently supported is GCC, CLANG (LLVM), and ICC (Intel). 1. Open `config.mk` and edit the `TAG` value according to the tool chain used. Currently supported is GCC, CLANG (LLVM), and ICC (Intel).
2. Open and adapt the compiler flags in `<include_<TOOLCHAIN>.mk`, e.g. in `include_ICC.mk` for the Intel tool chain. 2. Change `DATA_LAYOUT` and `DATA_TYPE` if desired in config.mk.
3. Build the binary calling `make`. 3. Open and adapt the compiler flags in `<include_<TOOLCHAIN>.mk`, e.g. in `include_ICC.mk` for the Intel tool chain.
4. Build the binary calling `make`.
You can clean intermediate build results with `make clean`, and all build results with `make distclean`. You can clean intermediate build results with `make clean`, and all build results with `make distclean`.
You have to call `make clean` before `make` if you changed the build settings. You have to call `make clean` before `make` if you changed the build settings.

7
config.mk Normal file
View File

@ -0,0 +1,7 @@
# Supported: GCC, CLANG, ICC
TAG ?= GCC
DATA_TYPE ?= SP#SP or DP
DATA_LAYOUT ?= SoA#AOS or SOA
#Feature options
OPTIONS += -DALIGNMENT=64

9
include_CLANG.mk
View File

@ -1,5 +1,4 @@
CC = cc CC = cc
CXX = cc
LINKER = $(CC) LINKER = $(CC)
ANSI_CFLAGS = -ansi ANSI_CFLAGS = -ansi
@ -7,11 +6,9 @@ ANSI_CFLAGS += -std=c99
ANSI_CFLAGS += -pedantic ANSI_CFLAGS += -pedantic
ANSI_CFLAGS += -Wextra ANSI_CFLAGS += -Wextra
CFLAGS = -Ofast $(ANSI_CFLAGS) -Xpreprocessor -fopenmp #-g CFLAGS = -Ofast $(ANSI_CFLAGS) -g #-Xpreprocessor -fopenmp -g
ASFLAGS = -masm=intel ASFLAGS = -masm=intel
CXXFLAGS = $(CFLAGS)
FCFLAGS =
LFLAGS = LFLAGS =
DEFINES = -D_GNU_SOURCE -DALIGNMENT=64 -DPRECISION=2 DEFINES = -D_GNU_SOURCE
INCLUDES = INCLUDES =
LIBS = -lomp LIBS = -lm #-lomp

12
include_GCC.mk
View File

@ -1,7 +1,5 @@
CC = gcc CC = gcc
CXX = g++ LINKER = $(CC)
FC = gfortran
LINKER = $(CXX)
ANSI_CFLAGS = -ansi ANSI_CFLAGS = -ansi
ANSI_CFLAGS += -std=c99 ANSI_CFLAGS += -std=c99
@ -9,13 +7,9 @@ ANSI_CFLAGS += -pedantic
ANSI_CFLAGS += -Wextra ANSI_CFLAGS += -Wextra
# CFLAGS = -O0 -g -std=c99 -fargument-noalias # CFLAGS = -O0 -g -std=c99 -fargument-noalias
CFLAGS = -O3 -march=znver1 -ffast-math -funroll-loops -fopenmp CFLAGS = -O3 -march=znver1 -ffast-math -funroll-loops # -fopenmp
CXXFLAGS = $(CFLAGS)
ASFLAGS = -masm=intel ASFLAGS = -masm=intel
FCFLAGS =
LFLAGS = LFLAGS =
DEFINES = -D_GNU_SOURCE DEFINES = -D_GNU_SOURCE
INCLUDES = INCLUDES =
LIBS = LIBS = -lm

7
include_ICC.mk
View File

@ -1,5 +1,4 @@
CC = icc CC = icc
CXX = icpc
LINKER = $(CC) LINKER = $(CC)
OPENMP = #-qopenmp OPENMP = #-qopenmp
@ -11,10 +10,8 @@ PROFILE = #-profile-functions -g -pg
#OPTS = -fast -no-vec $(PROFILE) #OPTS = -fast -no-vec $(PROFILE)
OPTS = -fast -xHost $(PROFILE) OPTS = -fast -xHost $(PROFILE)
CFLAGS = $(PROFILE) -restrict $(OPENMP) $(OPTS) CFLAGS = $(PROFILE) -restrict $(OPENMP) $(OPTS)
CXXFLAGS = $(CFLAGS)
ASFLAGS = -masm=intel ASFLAGS = -masm=intel
FCFLAGS =
LFLAGS = $(PROFILE) $(OPTS) $(OPENMP) LFLAGS = $(PROFILE) $(OPTS) $(OPENMP)
DEFINES = -D_GNU_SOURCE -DALIGNMENT=64 # -DLIKWID_PERFMON -DPRECISION=1 DEFINES = -D_GNU_SOURCE # -DALIGNMENT=64 -DLIKWID_PERFMON -DPRECISION=1
INCLUDES = #$(LIKWID_INC) INCLUDES = #$(LIKWID_INC)
LIBS = #$(LIKWID_LIB) -llikwid LIBS = -lm #$(LIKWID_LIB) -llikwid

16
src/atom.c
View File

@ -1,8 +1,10 @@
/* /*
* ======================================================================================= * =======================================================================================
* *
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de * Authors: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg * Rafael Ravedutti (rr), rafaelravedutti@gmail.com
*
* Copyright (c) 2021 RRZE, University Erlangen-Nuremberg
* *
* This file is part of MD-Bench. * This file is part of MD-Bench.
* *
@ -111,9 +113,9 @@ void createAtom(Atom *atom, Parameter *param)
growAtom(atom); growAtom(atom);
} }
atom->x[atom->Nlocal] = xtmp; atom_x(atom->Nlocal) = xtmp;
atom->y[atom->Nlocal] = ytmp; atom_y(atom->Nlocal) = ytmp;
atom->z[atom->Nlocal] = ztmp; atom_z(atom->Nlocal) = ztmp;
atom->vx[atom->Nlocal] = vxtmp; atom->vx[atom->Nlocal] = vxtmp;
atom->vy[atom->Nlocal] = vytmp; atom->vy[atom->Nlocal] = vytmp;
atom->vz[atom->Nlocal] = vztmp; atom->vz[atom->Nlocal] = vztmp;
@ -136,9 +138,13 @@ void growAtom(Atom *atom)
int nold = atom->Nmax; int nold = atom->Nmax;
atom->Nmax += DELTA; atom->Nmax += DELTA;
#ifdef AOS
atom->x = (MD_FLOAT*) reallocate(atom->x, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT) * 3, nold * sizeof(MD_FLOAT) * 3);
#else
atom->x = (MD_FLOAT*) reallocate(atom->x, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT)); atom->x = (MD_FLOAT*) reallocate(atom->x, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT));
atom->y = (MD_FLOAT*) reallocate(atom->y, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT)); atom->y = (MD_FLOAT*) reallocate(atom->y, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT));
atom->z = (MD_FLOAT*) reallocate(atom->z, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT)); atom->z = (MD_FLOAT*) reallocate(atom->z, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT));
#endif
atom->vx = (MD_FLOAT*) reallocate(atom->vx, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT)); atom->vx = (MD_FLOAT*) reallocate(atom->vx, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT));
atom->vy = (MD_FLOAT*) reallocate(atom->vy, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT)); atom->vy = (MD_FLOAT*) reallocate(atom->vy, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT));
atom->vz = (MD_FLOAT*) reallocate(atom->vz, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT)); atom->vz = (MD_FLOAT*) reallocate(atom->vz, ALIGNMENT, atom->Nmax * sizeof(MD_FLOAT), nold * sizeof(MD_FLOAT));

15
src/includes/atom.h
View File

@ -2,7 +2,7 @@
* ======================================================================================= * =======================================================================================
* *
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de * Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg * Copyright (c) 2021 RRZE, University Erlangen-Nuremberg
* *
* This file is part of MD-Bench. * This file is part of MD-Bench.
* *
@ -35,4 +35,17 @@ typedef struct {
extern void initAtom(Atom*); extern void initAtom(Atom*);
extern void createAtom(Atom*, Parameter*); extern void createAtom(Atom*, Parameter*);
extern void growAtom(Atom*); extern void growAtom(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]
#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]
#endif
#endif #endif

3
src/includes/parameter.h
View File

@ -23,9 +23,6 @@
#ifndef __PARAMETER_H_ #ifndef __PARAMETER_H_
#define __PARAMETER_H_ #define __PARAMETER_H_
#ifndef PRECISION
#define PRECISION 2
#endif
#if PRECISION == 1 #if PRECISION == 1
#define MD_FLOAT float #define MD_FLOAT float
#else #else

21
src/main.c
View File

@ -116,7 +116,6 @@ double reneighbour(
void initialIntegrate(Parameter *param, Atom *atom) void initialIntegrate(Parameter *param, Atom *atom)
{ {
MD_FLOAT* x = atom->x; MD_FLOAT* y = atom->y; MD_FLOAT* z = atom->z;
MD_FLOAT* fx = atom->fx; MD_FLOAT* fy = atom->fy; MD_FLOAT* fz = atom->fz; MD_FLOAT* fx = atom->fx; MD_FLOAT* fy = atom->fy; MD_FLOAT* fz = atom->fz;
MD_FLOAT* vx = atom->vx; MD_FLOAT* vy = atom->vy; MD_FLOAT* vz = atom->vz; MD_FLOAT* vx = atom->vx; MD_FLOAT* vy = atom->vy; MD_FLOAT* vz = atom->vz;
@ -124,9 +123,9 @@ void initialIntegrate(Parameter *param, Atom *atom)
vx[i] += param->dtforce * fx[i]; vx[i] += param->dtforce * fx[i];
vy[i] += param->dtforce * fy[i]; vy[i] += param->dtforce * fy[i];
vz[i] += param->dtforce * fz[i]; vz[i] += param->dtforce * fz[i];
x[i] += param->dt * vx[i]; atom_x(i) = atom_x(i) + param->dt * vx[i];
y[i] += param->dt * vy[i]; atom_y(i) = atom_y(i) + param->dt * vy[i];
z[i] += param->dt * vz[i]; atom_z(i) = atom_z(i) + param->dt * vz[i];
} }
} }
@ -149,7 +148,6 @@ double computeForce(Parameter *param, Atom *atom, Neighbor *neighbor)
MD_FLOAT cutforcesq = param->cutforce * param->cutforce; MD_FLOAT cutforcesq = param->cutforce * param->cutforce;
MD_FLOAT sigma6 = param->sigma6; MD_FLOAT sigma6 = param->sigma6;
MD_FLOAT epsilon = param->epsilon; MD_FLOAT epsilon = param->epsilon;
MD_FLOAT* x = atom->x; MD_FLOAT* y = atom->y; MD_FLOAT* z = atom->z;
MD_FLOAT* fx = atom->fx; MD_FLOAT* fy = atom->fy; MD_FLOAT* fz = atom->fz; MD_FLOAT* fx = atom->fx; MD_FLOAT* fy = atom->fy; MD_FLOAT* fz = atom->fz;
MD_FLOAT S, E; MD_FLOAT S, E;
@ -166,9 +164,9 @@ double computeForce(Parameter *param, Atom *atom, Neighbor *neighbor)
for(int i = 0; i < Nlocal; i++) { for(int i = 0; i < Nlocal; i++) {
neighs = &neighbor->neighbors[i * neighbor->maxneighs]; neighs = &neighbor->neighbors[i * neighbor->maxneighs];
int numneighs = neighbor->numneigh[i]; int numneighs = neighbor->numneigh[i];
MD_FLOAT xtmp = x[i]; MD_FLOAT xtmp = atom_x(i);
MD_FLOAT ytmp = y[i]; MD_FLOAT ytmp = atom_y(i);
MD_FLOAT ztmp = z[i]; MD_FLOAT ztmp = atom_z(i);
MD_FLOAT fix = 0; MD_FLOAT fix = 0;
MD_FLOAT fiy = 0; MD_FLOAT fiy = 0;
@ -176,9 +174,9 @@ double computeForce(Parameter *param, Atom *atom, Neighbor *neighbor)
for(int k = 0; k < numneighs; k++) { for(int k = 0; k < numneighs; k++) {
int j = neighs[k]; int j = neighs[k];
MD_FLOAT delx = xtmp - x[j]; MD_FLOAT delx = xtmp - atom_x(j);
MD_FLOAT dely = ytmp - y[j]; MD_FLOAT dely = ytmp - atom_y(j);
MD_FLOAT delz = ztmp - z[j]; MD_FLOAT delz = ztmp - atom_z(j);
MD_FLOAT rsq = delx * delx + dely * dely + delz * delz; MD_FLOAT rsq = delx * delx + dely * dely + delz * delz;
if(rsq < cutforcesq) { if(rsq < cutforcesq) {
@ -293,6 +291,7 @@ int main (int argc, char** argv)
computeThermo(-1, &param, &atom); computeThermo(-1, &param, &atom);
printf(HLINE); printf(HLINE);
printf("Data layout for positions: %s\n", POS_DATA_LAYOUT);
#if PRECISION == 1 #if PRECISION == 1
printf("Using single precision floating point.\n"); printf("Using single precision floating point.\n");
#else #else

22
src/neighbor.c
View File

@ -2,7 +2,7 @@
* ======================================================================================= * =======================================================================================
* *
* Author: Jan Eitzinger (je), jan.eitzinger@fau.de * Author: Jan Eitzinger (je), jan.eitzinger@fau.de
* Copyright (c) 2020 RRZE, University Erlangen-Nuremberg * Copyright (c) 2021 RRZE, University Erlangen-Nuremberg
* *
* This file is part of MD-Bench. * This file is part of MD-Bench.
* *
@ -184,9 +184,6 @@ void buildNeighbor(Atom *atom, Neighbor *neighbor)
/* bin local & ghost atoms */ /* bin local & ghost atoms */
binatoms(atom); binatoms(atom);
int resize = 1; int resize = 1;
MD_FLOAT* x = atom->x;
MD_FLOAT* y = atom->y;
MD_FLOAT* z = atom->z;
/* loop over each atom, storing neighbors */ /* loop over each atom, storing neighbors */
while(resize) { while(resize) {
@ -196,9 +193,9 @@ void buildNeighbor(Atom *atom, Neighbor *neighbor)
for(int i = 0; i < atom->Nlocal; i++) { for(int i = 0; i < atom->Nlocal; i++) {
int* neighptr = &(neighbor->neighbors[i * neighbor->maxneighs]); int* neighptr = &(neighbor->neighbors[i * neighbor->maxneighs]);
int n = 0; int n = 0;
MD_FLOAT xtmp = x[i]; MD_FLOAT xtmp = atom_x(i);
MD_FLOAT ytmp = y[i]; MD_FLOAT ytmp = atom_y(i);
MD_FLOAT ztmp = z[i]; MD_FLOAT ztmp = atom_z(i);
int ibin = coord2bin(xtmp, ytmp, ztmp); int ibin = coord2bin(xtmp, ytmp, ztmp);
for(int k = 0; k < nstencil; k++) { for(int k = 0; k < nstencil; k++) {
@ -212,9 +209,9 @@ void buildNeighbor(Atom *atom, Neighbor *neighbor)
continue; continue;
} }
MD_FLOAT delx = xtmp - x[j]; MD_FLOAT delx = xtmp - atom_x(j);
MD_FLOAT dely = ytmp - y[j]; MD_FLOAT dely = ytmp - atom_y(j);
MD_FLOAT delz = ztmp - z[j]; MD_FLOAT delz = ztmp - atom_z(j);
MD_FLOAT rsq = delx * delx + dely * dely + delz * delz; MD_FLOAT rsq = delx * delx + dely * dely + delz * delz;
if( rsq <= cutneighsq ) { if( rsq <= cutneighsq ) {
@ -309,9 +306,6 @@ int coord2bin(MD_FLOAT xin, MD_FLOAT yin, MD_FLOAT zin)
void binatoms(Atom *atom) void binatoms(Atom *atom)
{ {
int nall = atom->Nlocal + atom->Nghost; int nall = atom->Nlocal + atom->Nghost;
MD_FLOAT* x = atom->x;
MD_FLOAT* y = atom->y;
MD_FLOAT* z = atom->z;
int resize = 1; int resize = 1;
while(resize > 0) { while(resize > 0) {
@ -322,7 +316,7 @@ void binatoms(Atom *atom)
} }
for(int i = 0; i < nall; i++) { for(int i = 0; i < nall; i++) {
int ibin = coord2bin(x[i], y[i], z[i]); int ibin = coord2bin(atom_x(i), atom_y(i), atom_z(i));
if(bincount[ibin] < atoms_per_bin) { if(bincount[ibin] < atoms_per_bin) {
int ac = bincount[ibin]++; int ac = bincount[ibin]++;

101
src/pbc.c
View File

@ -48,17 +48,14 @@ void initPbc()
void updatePbc(Atom *atom, Parameter *param) void updatePbc(Atom *atom, Parameter *param)
{ {
int nlocal = atom->Nlocal; int nlocal = atom->Nlocal;
MD_FLOAT* x = atom->x;
MD_FLOAT* y = atom->y;
MD_FLOAT* z = atom->z;
MD_FLOAT xprd = param->xprd; MD_FLOAT xprd = param->xprd;
MD_FLOAT yprd = param->yprd; MD_FLOAT yprd = param->yprd;
MD_FLOAT zprd = param->zprd; MD_FLOAT zprd = param->zprd;
for(int i = 0; i < atom->Nghost; i++) { for(int i = 0; i < atom->Nghost; i++) {
x[nlocal + i] = x[BorderMap[i]] + PBCx[i] * xprd; atom_x(nlocal + i) = atom_x(BorderMap[i]) + PBCx[i] * xprd;
y[nlocal + i] = y[BorderMap[i]] + PBCy[i] * yprd; atom_y(nlocal + i) = atom_y(BorderMap[i]) + PBCy[i] * yprd;
z[nlocal + i] = z[BorderMap[i]] + PBCz[i] * zprd; atom_z(nlocal + i) = atom_z(BorderMap[i]) + PBCz[i] * zprd;
} }
} }
@ -66,31 +63,28 @@ void updatePbc(Atom *atom, Parameter *param)
* to periodic boundary conditions */ * to periodic boundary conditions */
void updateAtomsPbc(Atom *atom, Parameter *param) void updateAtomsPbc(Atom *atom, Parameter *param)
{ {
MD_FLOAT* x = atom->x;
MD_FLOAT* y = atom->y;
MD_FLOAT* z = atom->z;
MD_FLOAT xprd = param->xprd; MD_FLOAT xprd = param->xprd;
MD_FLOAT yprd = param->yprd; MD_FLOAT yprd = param->yprd;
MD_FLOAT zprd = param->zprd; MD_FLOAT zprd = param->zprd;
for(int i = 0; i < atom->Nlocal; i++) { for(int i = 0; i < atom->Nlocal; i++) {
if(x[i] < 0.0) { if(atom_x(i) < 0.0) {
x[i] += xprd; atom_x(i) += xprd;
} else if(x[i] >= xprd) { } else if(atom_x(i) >= xprd) {
x[i] -= xprd; atom_x(i) -= xprd;
} }
if(y[i] < 0.0) { if(atom_y(i) < 0.0) {
y[i] += yprd; atom_y(i) += yprd;
} else if(y[i] >= yprd) { } else if(atom_y(i) >= yprd) {
y[i] -= yprd; atom_y(i) -= yprd;
} }
if(z[i] < 0.0) { if(atom_z(i) < 0.0) {
z[i] += zprd; atom_z(i) += zprd;
} else if(z[i] >= zprd) { } else if(atom_z(i) >= zprd) {
z[i] -= zprd; atom_z(i) -= zprd;
} }
} }
} }
@ -99,10 +93,14 @@ 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) Nghost++; BorderMap[Nghost] = i; PBCx[Nghost] = dx; PBCy[Nghost] = dy; PBCz[Nghost] = dz; #define ADDGHOST(dx,dy,dz) \
Nghost++; \
BorderMap[Nghost] = i; \
PBCx[Nghost] = dx; \
PBCy[Nghost] = dy; \
PBCz[Nghost] = dz;
void setupPbc(Atom *atom, Parameter *param) void setupPbc(Atom *atom, Parameter *param)
{ {
MD_FLOAT* x = atom->x; MD_FLOAT* y = atom->y; MD_FLOAT* z = atom->z;
MD_FLOAT xprd = param->xprd; MD_FLOAT xprd = param->xprd;
MD_FLOAT yprd = param->yprd; MD_FLOAT yprd = param->yprd;
MD_FLOAT zprd = param->zprd; MD_FLOAT zprd = param->zprd;
@ -113,42 +111,45 @@ void setupPbc(Atom *atom, Parameter *param)
if (atom->Nlocal + Nghost + 7 >= atom->Nmax) { if (atom->Nlocal + Nghost + 7 >= atom->Nmax) {
growAtom(atom); growAtom(atom);
x = atom->x; y = atom->y; z = atom->z;
} }
if (Nghost + 7 >= NmaxGhost) { if (Nghost + 7 >= NmaxGhost) {
growPbc(); growPbc();
} }
MD_FLOAT x = atom_x(i);
MD_FLOAT y = atom_y(i);
MD_FLOAT z = atom_z(i);
/* Setup ghost atoms */ /* Setup ghost atoms */
/* 6 planes */ /* 6 planes */
if (x[i] < Cutneigh) { ADDGHOST(+1,0,0); } if (x < Cutneigh) { ADDGHOST(+1,0,0); }
if (x[i] >= (xprd-Cutneigh)) { ADDGHOST(-1,0,0); } if (x >= (xprd-Cutneigh)) { ADDGHOST(-1,0,0); }
if (y[i] < Cutneigh) { ADDGHOST(0,+1,0); } if (y < Cutneigh) { ADDGHOST(0,+1,0); }
if (y[i] >= (yprd-Cutneigh)) { ADDGHOST(0,-1,0); } if (y >= (yprd-Cutneigh)) { ADDGHOST(0,-1,0); }
if (z[i] < Cutneigh) { ADDGHOST(0,0,+1); } if (z < Cutneigh) { ADDGHOST(0,0,+1); }
if (z[i] >= (zprd-Cutneigh)) { ADDGHOST(0,0,-1); } if (z >= (zprd-Cutneigh)) { ADDGHOST(0,0,-1); }
/* 8 corners */ /* 8 corners */
if (x[i] < Cutneigh && y[i] < Cutneigh && z[i] < Cutneigh) { ADDGHOST(+1,+1,+1); } if (x < Cutneigh && y < Cutneigh && z < Cutneigh) { ADDGHOST(+1,+1,+1); }
if (x[i] < Cutneigh && y[i] >= (yprd-Cutneigh) && z[i] < Cutneigh) { ADDGHOST(+1,-1,+1); } if (x < Cutneigh && y >= (yprd-Cutneigh) && z < Cutneigh) { ADDGHOST(+1,-1,+1); }
if (x[i] < Cutneigh && y[i] >= Cutneigh && z[i] >= (zprd-Cutneigh)) { ADDGHOST(+1,+1,-1); } if (x < Cutneigh && y >= Cutneigh && z >= (zprd-Cutneigh)) { ADDGHOST(+1,+1,-1); }
if (x[i] < Cutneigh && y[i] >= (yprd-Cutneigh) && z[i] >= (zprd-Cutneigh)) { ADDGHOST(+1,-1,-1); } if (x < Cutneigh && y >= (yprd-Cutneigh) && z >= (zprd-Cutneigh)) { ADDGHOST(+1,-1,-1); }
if (x[i] >= (xprd-Cutneigh) && y[i] < Cutneigh && z[i] < Cutneigh) { ADDGHOST(-1,+1,+1); } if (x >= (xprd-Cutneigh) && y < Cutneigh && z < Cutneigh) { ADDGHOST(-1,+1,+1); }
if (x[i] >= (xprd-Cutneigh) && y[i] >= (yprd-Cutneigh) && z[i] < Cutneigh) { ADDGHOST(-1,-1,+1); } if (x >= (xprd-Cutneigh) && y >= (yprd-Cutneigh) && z < Cutneigh) { ADDGHOST(-1,-1,+1); }
if (x[i] >= (xprd-Cutneigh) && y[i] < Cutneigh && z[i] >= (zprd-Cutneigh)) { ADDGHOST(-1,+1,-1); } if (x >= (xprd-Cutneigh) && y < Cutneigh && z >= (zprd-Cutneigh)) { ADDGHOST(-1,+1,-1); }
if (x[i] >= (xprd-Cutneigh) && y[i] >= (yprd-Cutneigh) && z[i] >= (zprd-Cutneigh)) { ADDGHOST(-1,-1,-1); } if (x >= (xprd-Cutneigh) && y >= (yprd-Cutneigh) && z >= (zprd-Cutneigh)) { ADDGHOST(-1,-1,-1); }
/* 12 edges */ /* 12 edges */
if (x[i] < Cutneigh && z[i] < Cutneigh) { ADDGHOST(+1,0,+1); } if (x < Cutneigh && z < Cutneigh) { ADDGHOST(+1,0,+1); }
if (x[i] < Cutneigh && z[i] >= (zprd-Cutneigh)) { ADDGHOST(+1,0,-1); } if (x < Cutneigh && z >= (zprd-Cutneigh)) { ADDGHOST(+1,0,-1); }
if (x[i] >= (xprd-Cutneigh) && z[i] < Cutneigh) { ADDGHOST(-1,0,+1); } if (x >= (xprd-Cutneigh) && z < Cutneigh) { ADDGHOST(-1,0,+1); }
if (x[i] >= (xprd-Cutneigh) && z[i] >= (zprd-Cutneigh)) { ADDGHOST(-1,0,-1); } if (x >= (xprd-Cutneigh) && z >= (zprd-Cutneigh)) { ADDGHOST(-1,0,-1); }
if (y[i] < Cutneigh && z[i] < Cutneigh) { ADDGHOST(0,+1,+1); } if (y < Cutneigh && z < Cutneigh) { ADDGHOST(0,+1,+1); }
if (y[i] < Cutneigh && z[i] >= (zprd-Cutneigh)) { ADDGHOST(0,+1,-1); } if (y < Cutneigh && z >= (zprd-Cutneigh)) { ADDGHOST(0,+1,-1); }
if (y[i] >= (yprd-Cutneigh) && z[i] < Cutneigh) { ADDGHOST(0,-1,+1); } if (y >= (yprd-Cutneigh) && z < Cutneigh) { ADDGHOST(0,-1,+1); }
if (y[i] >= (yprd-Cutneigh) && z[i] >= (zprd-Cutneigh)) { ADDGHOST(0,-1,-1); } if (y >= (yprd-Cutneigh) && z >= (zprd-Cutneigh)) { ADDGHOST(0,-1,-1); }
if (y[i] < Cutneigh && x[i] < Cutneigh) { ADDGHOST(+1,+1,0); } if (y < Cutneigh && x < Cutneigh) { ADDGHOST(+1,+1,0); }
if (y[i] < Cutneigh && x[i] >= (xprd-Cutneigh)) { ADDGHOST(-1,+1,0); } if (y < Cutneigh && x >= (xprd-Cutneigh)) { ADDGHOST(-1,+1,0); }
if (y[i] >= (yprd-Cutneigh) && x[i] < Cutneigh) { ADDGHOST(+1,-1,0); } if (y >= (yprd-Cutneigh) && x < Cutneigh) { ADDGHOST(+1,-1,0); }
if (y[i] >= (yprd-Cutneigh) && x[i] >= (xprd-Cutneigh)) { ADDGHOST(-1,-1,0); } if (y >= (yprd-Cutneigh) && x >= (xprd-Cutneigh)) { ADDGHOST(-1,-1,0); }
} }
// increase by one to make it the ghost atom count // increase by one to make it the ghost atom count
atom->Nghost = Nghost + 1; atom->Nghost = Nghost + 1;