getting there

This commit is contained in:
Tom Swinburne
2020-04-13 19:28:43 +02:00
parent 39799c62fc
commit 430f2ae6aa
9 changed files with 376 additions and 2212 deletions

View File

@ -1490,24 +1490,8 @@ void lammps_scatter_atoms_subset(void *ptr, char *name,
#endif
/* ----------------------------------------------------------------------
Contributing author: Thomas Swinburne (CNRS & CINaM, Marseille, France)
gather the named per atom fix and return it in user-allocated data
(extract_fix was found to give errors when running lammps in parallel)
data will be ordered by atom ID
requirement for consecutive atom IDs (1 to N)
id = fix ID
count: number of entries per atom
Fills 1d data, which must be pre-allocated to length of count * Natoms, where Natoms is as queried by get_natoms()
method:
alloc and zero count*Natom length vector
loop over nlocal to fill vector with my values
Allreduce to sum vector into data across all procs
------------------------------------------------------------------------- */
#if defined(LAMMPS_BIGBIG)
void lammps_gather_peratom_fix(void *ptr, char * /*id */, int /*count*/, void * /*data*/)
void lammps_gather_peratom_fix(void *ptr, char * /*id */, int /*type*/, int /*count*/, void * /*data*/)
{
LAMMPS *lmp = (LAMMPS *) ptr;
@ -1516,7 +1500,7 @@ void lammps_gather_peratom_fix(void *ptr, char * /*id */, int /*count*/, void *
END_CAPTURE
}
#else
void lammps_gather_peratom_fix(void *ptr, char *id, int count, void *data)
void lammps_gather_peratom_fix(void *ptr, char *id, int type, int count, void *data)
{
LAMMPS *lmp = (LAMMPS *) ptr;
@ -1526,7 +1510,7 @@ void lammps_gather_peratom_fix(void *ptr, char *id, int count, void *data)
int ifix = lmp->modify->find_fix(id);
if (ifix < 0) {
lmp->error->warning(FLERR,"lammps_gather_pertatom_fix: unknown fix id");
lmp->error->warning(FLERR,"lammps_gather_peratom_fix: unknown fix id");
return;
}
@ -1547,31 +1531,57 @@ void lammps_gather_peratom_fix(void *ptr, char *id, int count, void *data)
// copy = Natom length vector of per-atom values
// use atom ID to insert each atom's values into copy
// MPI_Allreduce with MPI_SUM to merge into data, ordered by atom ID
if (type==0) {
int *vector = NULL;
int **array = NULL;
if (count == 1) vector = (int *) fix->vector_atom; //vptr;
else array = (int **) fix->array_atom; //vptr;
double *vector = NULL;
double **array = NULL;
if (count == 1) vector = (double *) fix->vector_atom; //vptr;
else array = (double **) fix->array_atom; //vptr;
int *copy;
lmp->memory->create(copy,count*natoms,"lib/gather:copy");
for (i = 0; i < count*natoms; i++) copy[i] = 0;
double *copy;
lmp->memory->create(copy,count*natoms,"lib/gather:copy");
for (i = 0; i < count*natoms; i++) copy[i] = 0.0;
tagint *tag = lmp->atom->tag;
int nlocal = lmp->atom->nlocal;
tagint *tag = lmp->atom->tag;
int nlocal = lmp->atom->nlocal;
if (count == 1) {
for (i = 0; i < nlocal; i++)
copy[tag[i]-1] = vector[i];
} else {
for (i = 0; i < nlocal; i++) {
offset = count*(tag[i]-1);
for (j = 0; j < count; j++)
copy[offset++] = array[i][j];
if (count == 1) {
for (i = 0; i < nlocal; i++)
copy[tag[i]-1] = vector[i];
} else {
for (i = 0; i < nlocal; i++) {
offset = count*(tag[i]-1);
for (j = 0; j < count; j++)
copy[offset++] = array[i][j];
}
}
MPI_Allreduce(copy,data,count*natoms,MPI_INT,MPI_SUM,lmp->world);
lmp->memory->destroy(copy);
} else {
double *vector = NULL;
double **array = NULL;
if (count == 1) vector = (double *) fix->vector_atom; //vptr;
else array = (double **) fix->array_atom; //vptr;
double *copy;
lmp->memory->create(copy,count*natoms,"lib/gather:copy");
for (i = 0; i < count*natoms; i++) copy[i] = 0.0;
tagint *tag = lmp->atom->tag;
int nlocal = lmp->atom->nlocal;
if (count == 1) {
for (i = 0; i < nlocal; i++)
copy[tag[i]-1] = vector[i];
} else {
for (i = 0; i < nlocal; i++) {
offset = count*(tag[i]-1);
for (j = 0; j < count; j++)
copy[offset++] = array[i][j];
}
}
MPI_Allreduce(copy,data,count*natoms,MPI_DOUBLE,MPI_SUM,lmp->world);
lmp->memory->destroy(copy);
}
MPI_Allreduce(copy,data,count*natoms,MPI_DOUBLE,MPI_SUM,lmp->world);
lmp->memory->destroy(copy);
}
END_CAPTURE
}
@ -1595,8 +1605,8 @@ void lammps_gather_peratom_fix(void *ptr, char *id, int count, void *data)
Allreduce to sum vector into data across all procs
------------------------------------------------------------------------- */
#if defined(LAMMPS_BIGBIG)
void lammps_gather_peratom_fix_subset(void *ptr, char * /*id */, int /*count*/,
int /*ndata*/, int * /*ids*/, void * /*data*/)
void lammps_gather_peratom_fix_subset(void *ptr, char * /*id */, int /*type*/,
int /*count*/, int /*ndata*/, int * /*ids*/, void * /*data*/)
{
LAMMPS *lmp = (LAMMPS *) ptr;
@ -1605,7 +1615,7 @@ void lammps_gather_peratom_fix_subset(void *ptr, char * /*id */, int /*count*/,
END_CAPTURE
}
#else
void lammps_gather_peratom_fix_subset(void *ptr, char *id, int count,
void lammps_gather_peratom_fix_subset(void *ptr, char *id, int type, int count,
int ndata, int *ids, void *data)
{
LAMMPS *lmp = (LAMMPS *) ptr;
@ -1617,7 +1627,7 @@ void lammps_gather_peratom_fix_subset(void *ptr, char *id, int count,
int ifix = lmp->modify->find_fix(id);
if (ifix < 0) {
lmp->error->warning(FLERR,"lammps_gather_pertatom_fix: unknown fix id");
lmp->error->warning(FLERR,"lammps_gather_peratom_fix_subset: unknown fix id");
return;
}
@ -1631,49 +1641,330 @@ void lammps_gather_peratom_fix_subset(void *ptr, char *id, int count,
if (lmp->atom->natoms > MAXSMALLINT) flag = 1;
if (flag) {
if (lmp->comm->me == 0)
lmp->error->warning(FLERR,"Library error in lammps_gather_peratom_fix");
lmp->error->warning(FLERR,"Library error in lammps_gather_peratom_fix_subset");
return;
}
// copy = Natom length vector of per-atom values
// use atom ID to insert each atom's values into copy
// MPI_Allreduce with MPI_SUM to merge into data, ordered by atom ID
if (type==0) {
int *vector = NULL;
int **array = NULL;
if (count == 1) vector = (int *) fix->vector_atom; //vptr;
else array = (int **) fix->array_atom; //vptr;
double *vector = NULL;
double **array = NULL;
if (count == 1) vector = (double *) fix->vector_atom; //vptr;
else array = (double **) fix->array_atom; //vptr;
int *copy;
lmp->memory->create(copy,count*ndata,"lib/gather:copy");
for (i = 0; i < count*ndata; i++) copy[i] = 0;
double *copy;
lmp->memory->create(copy,count*natoms,"lib/gather:copy");
for (i = 0; i < count*natoms; i++) copy[i] = 0.0;
tagint *tag = lmp->atom->tag;
int nlocal = lmp->atom->nlocal;
tagint *tag = lmp->atom->tag;
int nlocal = lmp->atom->nlocal;
if (count == 1) {
for (i = 0; i < ndata; i++) {
aid = ids[i];
if ((m = lmp->atom->map(aid)) >= 0 && m < nlocal)
copy[i] = vector[m];
}
} else {
for (i = 0; i < ndata; i++) {
aid = ids[i];
if ((m = lmp->atom->map(aid)) >= 0 && m < nlocal) {
offset = count*i;
for (j = 0; j < count; j++)
copy[offset++] = array[m][j];
if (count == 1) {
for (i = 0; i < ndata; i++) {
aid = ids[i];
if ((m = lmp->atom->map(aid)) >= 0 && m < nlocal)
copy[i] = vector[m];
}
} else {
for (i = 0; i < ndata; i++) {
aid = ids[i];
if ((m = lmp->atom->map(aid)) >= 0 && m < nlocal) {
offset = count*i;
for (j = 0; j < count; j++)
copy[offset++] = array[m][j];
}
}
}
MPI_Allreduce(copy,data,count*ndata,MPI_INT,MPI_SUM,lmp->world);
lmp->memory->destroy(copy);
} else {
double *vector = NULL;
double **array = NULL;
if (count == 1) vector = (double *) fix->vector_atom; //vptr;
else array = (double **) fix->array_atom; //vptr;
double *copy;
lmp->memory->create(copy,count*ndata,"lib/gather:copy");
for (i = 0; i < count*ndata; i++) copy[i] = 0.0;
tagint *tag = lmp->atom->tag;
int nlocal = lmp->atom->nlocal;
if (count == 1) {
for (i = 0; i < ndata; i++) {
aid = ids[i];
if ((m = lmp->atom->map(aid)) >= 0 && m < nlocal)
copy[i] = vector[m];
}
} else {
for (i = 0; i < ndata; i++) {
aid = ids[i];
if ((m = lmp->atom->map(aid)) >= 0 && m < nlocal) {
offset = count*i;
for (j = 0; j < count; j++)
copy[offset++] = array[m][j];
}
}
}
MPI_Allreduce(copy,data,count*ndata,MPI_DOUBLE,MPI_SUM,lmp->world);
lmp->memory->destroy(copy);
}
MPI_Allreduce(copy,data,count*natoms,MPI_DOUBLE,MPI_SUM,lmp->world);
lmp->memory->destroy(copy);
}
END_CAPTURE
}
#endif
/* ----------------------------------------------------------------------
scatter the named atom-based entity in data to all atoms
data is ordered by atom ID
requirement for consecutive atom IDs (1 to N)
see scatter_atoms_subset() to scatter data for some (or all) atoms, unordered
name = desired quantity, e.g. x or charge
type = 0 for integer values, 1 for double values
count = # of per-atom values, e.g. 1 for type or charge, 3 for x or f
use count = 3 with "image" for xyz to be packed into single image flag
data = atom-based values in 1d data, ordered by count, then by atom ID
e.g. x[0][0],x[0][1],x[0][2],x[1][0],x[1][1],x[1][2],x[2][0],...
data must be correct length = count*Natoms, as queried by get_natoms()
method:
loop over Natoms, if I own atom ID, set its values from data
------------------------------------------------------------------------- */
#if defined(LAMMPS_BIGBIG)
void lammps_scatter_peratom_fix(void *ptr, char * /*id */, int /*type*/, int /*count*/, void * /*data*/)
{
LAMMPS *lmp = (LAMMPS *) ptr;
BEGIN_CAPTURE
{
lmp->error->all(FLERR,"Library function lammps_scatter_peratom_fix() "
"is not compatible with -DLAMMPS_BIGBIG");
}
END_CAPTURE
}
#else
void lammps_scatter_peratom_fix(void *ptr, char *id, int type, int count, void *data)
{
LAMMPS *lmp = (LAMMPS *) ptr;
BEGIN_CAPTURE
{
int i,j,m,offset;
int ifix = lmp->modify->find_fix(id);
if (ifix < 0) {
lmp->error->warning(FLERR,"lammps_scatter_peratom_fix: unknown fix id");
return;
}
// error if tags are not defined or not consecutive or no atom map
// NOTE: test that name = image or ids is not a 64-bit int in code?
int flag = 0;
if (lmp->atom->tag_enable == 0 || lmp->atom->tag_consecutive() == 0)
flag = 1;
if (lmp->atom->natoms > MAXSMALLINT) flag = 1;
if (lmp->atom->map_style == 0) flag = 1;
if (flag) {
if (lmp->comm->me == 0)
lmp->error->warning(FLERR,"Library error in lammps_scatter_peratom_fix");
return;
}
Fix *fix = lmp->modify->fix[ifix];
int natoms = static_cast<int> (lmp->atom->natoms);
if (type == 0) {
int *vector = NULL;
int **array = NULL;
if (count == 1) vector = (int *) fix->vector_atom;
else array = (int **) fix->array_atom;
int *dptr = (int *) data;
if (count == 1) {
for (i = 0; i < natoms; i++)
if ((m = lmp->atom->map(i+1)) >= 0)
vector[m] = dptr[i];
} else {
for (i = 0; i < natoms; i++)
if ((m = lmp->atom->map(i+1)) >= 0) {
offset = count*i;
for (j = 0; j < count; j++)
array[m][j] = dptr[offset++];
}
}
} else {
double *vector = NULL;
double **array = NULL;
if (count == 1) vector = (double *) fix->vector_atom;
else array = (double **) fix->array_atom;
double *dptr = (double *) data;
if (count == 1) {
for (i = 0; i < natoms; i++)
if ((m = lmp->atom->map(i+1)) >= 0)
vector[m] = dptr[i];
} else {
for (i = 0; i < natoms; i++) {
if ((m = lmp->atom->map(i+1)) >= 0) {
offset = count*i;
for (j = 0; j < count; j++)
array[m][j] = dptr[offset++];
}
}
}
}
}
END_CAPTURE
}
#endif
/* ----------------------------------------------------------------------
scatter the named atom-based entity in data to a subset of atoms
data is ordered by provided atom IDs
no requirement for consecutive atom IDs (1 to N)
see scatter_atoms() to scatter data for all atoms, ordered by consecutive IDs
name = desired quantity, e.g. x or charge
type = 0 for integer values, 1 for double values
count = # of per-atom values, e.g. 1 for type or charge, 3 for x or f
use count = 3 with "image" for xyz to be packed into single image flag
ndata = # of atoms in ids and data (could be all atoms)
ids = list of Ndata atom IDs to scatter data to
data = atom-based values in 1d data, ordered by count, then by atom ID
e.g. x[0][0],x[0][1],x[0][2],x[1][0],x[1][1],x[1][2],x[2][0],...
data must be correct length = count*Ndata
method:
loop over Ndata, if I own atom ID, set its values from data
------------------------------------------------------------------------- */
#if defined(LAMMPS_BIGBIG)
void lammps_scatter_peratom_fix_subset(void *ptr, char * /*name */,
int /*type*/, int /*count*/,
int /*ndata*/, int * /*ids*/, void * /*data*/)
{
LAMMPS *lmp = (LAMMPS *) ptr;
BEGIN_CAPTURE
{
lmp->error->all(FLERR,"Library function lammps_scatter_atoms_subset() "
"is not compatible with -DLAMMPS_BIGBIG");
}
END_CAPTURE
}
#else
void lammps_scatter_peratom_fix_subset(void *ptr, char *name,
int type, int count,
int ndata, int *ids, void *data)
{
LAMMPS *lmp = (LAMMPS *) ptr;
BEGIN_CAPTURE
{
int i,j,m,offset;
tagint id;
// error if tags are not defined or no atom map
// NOTE: test that name = image or ids is not a 64-bit int in code?
int flag = 0;
if (lmp->atom->tag_enable == 0) flag = 1;
if (lmp->atom->natoms > MAXSMALLINT) flag = 1;
if (lmp->atom->map_style == 0) flag = 1;
if (flag) {
if (lmp->comm->me == 0)
lmp->error->warning(FLERR,"Library error in lammps_scatter_atoms_subset");
return;
}
void *vptr = lmp->atom->extract(name);
if(vptr == NULL) {
lmp->error->warning(FLERR,
"lammps_scatter_atoms_subset: unknown property name");
return;
}
// copy = Natom length vector of per-atom values
// use atom ID to insert each atom's values into copy
// MPI_Allreduce with MPI_SUM to merge into data, ordered by atom ID
if (type == 0) {
int *vector = NULL;
int **array = NULL;
const int imgpack = (count == 3) && (strcmp(name,"image") == 0);
if ((count == 1) || imgpack) vector = (int *) vptr;
else array = (int **) vptr;
int *dptr = (int *) data;
if (count == 1) {
for (i = 0; i < ndata; i++) {
id = ids[i];
if ((m = lmp->atom->map(id)) >= 0)
vector[m] = dptr[i];
}
} else if (imgpack) {
for (i = 0; i < ndata; i++) {
id = ids[i];
if ((m = lmp->atom->map(id)) >= 0) {
offset = count*i;
int image = dptr[offset++] + IMGMAX;
image += (dptr[offset++] + IMGMAX) << IMGBITS;
image += (dptr[offset++] + IMGMAX) << IMG2BITS;
vector[m] = image;
}
}
} else {
for (i = 0; i < ndata; i++) {
id = ids[i];
if ((m = lmp->atom->map(id)) >= 0) {
offset = count*i;
for (j = 0; j < count; j++)
array[m][j] = dptr[offset++];
}
}
}
} else {
double *vector = NULL;
double **array = NULL;
if (count == 1) vector = (double *) vptr;
else array = (double **) vptr;
double *dptr = (double *) data;
if (count == 1) {
for (i = 0; i < ndata; i++) {
id = ids[i];
if ((m = lmp->atom->map(id)) >= 0)
vector[m] = dptr[i];
}
} else {
for (i = 0; i < ndata; i++) {
id = ids[i];
if ((m = lmp->atom->map(id)) >= 0) {
offset = count*i;
for (j = 0; j < count; j++)
array[m][j] = dptr[offset++];
}
}
}
}
}
END_CAPTURE
}
#endif
/* ----------------------------------------------------------------------
create N atoms and assign them to procs based on coords
id = atom IDs (optional, NULL will generate 1 to N)