refactor compute global/atom to use array of struct

This commit is contained in:
Axel Kohlmeyer
2022-10-05 06:33:08 -04:00
parent 874e1522e7
commit 2bc0825e4a
4 changed files with 203 additions and 224 deletions

View File

@ -153,7 +153,7 @@ void ComputeChunkSpreadAtom::init()
{ {
init_chunk(); init_chunk();
// set indices of all computes,fixes,variables // store references of all computes and fixes
for (auto &val : values) { for (auto &val : values) {
if (val.which == ArgInfo::COMPUTE) { if (val.which == ArgInfo::COMPUTE) {

View File

@ -26,19 +26,12 @@
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
enum{VECTOR,ARRAY};
#define BIG 1.0e20
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) : ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) :
Compute(lmp, narg, arg), Compute(lmp, narg, arg), indices(nullptr), varatom(nullptr), vecglobal(nullptr)
idref(nullptr), which(nullptr), argindex(nullptr), value2index(nullptr), ids(nullptr),
indices(nullptr), varatom(nullptr), vecglobal(nullptr)
{ {
if (narg < 5) error->all(FLERR,"Illegal compute global/atom command"); if (narg < 5) utils::missing_cmd_args(FLERR,"compute global/atom", error);
peratom_flag = 1; peratom_flag = 1;
@ -47,13 +40,13 @@ ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) :
int iarg = 3; int iarg = 3;
ArgInfo argi(arg[iarg]); ArgInfo argi(arg[iarg]);
whichref = argi.get_type(); reference.which = argi.get_type();
indexref = argi.get_index1(); reference.argindex = argi.get_index1();
idref = argi.copy_name(); reference.id = argi.get_name();
if ((whichref == ArgInfo::UNKNOWN) || (whichref == ArgInfo::NONE) if ((reference.which == ArgInfo::UNKNOWN) || (reference.which == ArgInfo::NONE)
|| (argi.get_dim() > 1)) || (argi.get_dim() > 1))
error->all(FLERR,"Illegal compute global/atom command"); error->all(FLERR,"Illegal compute global/atom index property: {}", arg[iarg]);
iarg++; iarg++;
@ -68,24 +61,21 @@ ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) :
// parse values until one isn't recognized // parse values until one isn't recognized
which = new int[nargnew]; values.clear();
argindex = new int[nargnew];
ids = new char*[nargnew];
value2index = new int[nargnew];
nvalues = 0;
for (iarg = 0; iarg < nargnew; iarg++) { for (iarg = 0; iarg < nargnew; iarg++) {
ArgInfo argi2(arg[iarg]); ArgInfo argi2(arg[iarg]);
which[nvalues] = argi2.get_type(); value_t val;
argindex[nvalues] = argi2.get_index1(); val.which = argi2.get_type();
ids[nvalues] = argi2.copy_name(); val.argindex = argi2.get_index1();
val.id = argi2.get_name();
if ((which[nvalues] == ArgInfo::UNKNOWN) || (which[nvalues] == ArgInfo::NONE) if ((val.which == ArgInfo::UNKNOWN) || (val.which == ArgInfo::NONE)
|| (argi2.get_dim() > 1)) || (argi2.get_dim() > 1))
error->all(FLERR,"Illegal compute slice command"); error->all(FLERR,"Illegal compute global/atom global property: {}", arg[iarg]);
nvalues++; values.push_back(val);
} }
// if wildcard expansion occurred, free earg memory from expand_args() // if wildcard expansion occurred, free earg memory from expand_args()
@ -95,94 +85,97 @@ ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) :
memory->sfree(earg); memory->sfree(earg);
} }
// setup and error check both index arg and values // setup and error check for both, index arg and values
if (whichref == ArgInfo::COMPUTE) { if (reference.which == ArgInfo::COMPUTE) {
int icompute = modify->find_compute(idref); reference.val.c = modify->get_compute_by_id(reference.id);
if (icompute < 0) if (!reference.val.c)
error->all(FLERR,"Compute ID for compute global/atom does not exist"); error->all(FLERR,"Compute ID {} for compute global/atom index", reference.id);
if (!modify->compute[icompute]->peratom_flag) if (!reference.val.c->peratom_flag)
error->all(FLERR,"Compute global/atom compute does not " error->all(FLERR,"Compute global/atom compute {} does not calculate a per-atom "
"calculate a per-atom vector or array"); "vector or array", reference.id);
if (indexref == 0 && if ((reference.argindex == 0) && (reference.val.c->size_peratom_cols != 0))
modify->compute[icompute]->size_peratom_cols != 0) error->all(FLERR,"Compute global/atom compute {} does not calculate a per-atom "
error->all(FLERR,"Compute global/atom compute does not " "vector", reference.id);
"calculate a per-atom vector"); if (reference.argindex && (reference.val.c->size_peratom_cols == 0))
if (indexref && modify->compute[icompute]->size_peratom_cols == 0) error->all(FLERR,"Compute global/atom compute does not calculate a per-atom "
error->all(FLERR,"Compute global/atom compute does not " "array", reference.id);
"calculate a per-atom array"); if (reference.argindex && (reference.argindex > reference.val.c->size_peratom_cols))
if (indexref && indexref > modify->compute[icompute]->size_peratom_cols) error->all(FLERR, "Compute global/atom compute array {} is accessed out-of-range",
error->all(FLERR, reference.id);
"Compute global/atom compute array is accessed out-of-range");
} else if (whichref == ArgInfo::FIX) { } else if (reference.which == ArgInfo::FIX) {
auto ifix = modify->get_fix_by_id(idref); reference.val.f =modify->get_fix_by_id(reference.id);
if (!ifix) if (!reference.val.f)
error->all(FLERR,"Fix ID {} for compute global/atom does not exist", idref); error->all(FLERR,"Fix ID {} for compute global/atom does not exist", reference.id);
if (!ifix->peratom_flag) if (!reference.val.f->peratom_flag)
error->all(FLERR,"Compute global/atom fix {} does not calculate a per-atom vector or array", idref); error->all(FLERR,"Compute global/atom fix {} does not calculate a per-atom vector "
if (indexref == 0 && (ifix->size_peratom_cols != 0)) "or array", reference.id);
error->all(FLERR,"Compute global/atom fix {} does not calculate a per-atom vector", idref); if (reference.argindex == 0 && (reference.val.f->size_peratom_cols != 0))
if (indexref && (ifix->size_peratom_cols == 0)) error->all(FLERR,"Compute global/atom fix {} does not calculate a per-atom vector",
error->all(FLERR,"Compute global/atom fix {} does not calculate a per-atom array", idref); reference.id);
if (indexref && indexref > ifix->size_peratom_cols) if (reference.argindex && (reference.val.f->size_peratom_cols == 0))
error->all(FLERR, "Compute global/atom fix {} array is accessed out-of-range", idref); error->all(FLERR,"Compute global/atom fix {} does not calculate a per-atom array",
reference.id);
if (reference.argindex && (reference.argindex > reference.val.f->size_peratom_cols))
error->all(FLERR, "Compute global/atom fix {} array is accessed out-of-range", reference.id);
} else if (whichref == ArgInfo::VARIABLE) { } else if (reference.which == ArgInfo::VARIABLE) {
int ivariable = input->variable->find(idref); reference.val.v = input->variable->find(reference.id.c_str());
if (ivariable < 0) if (reference.val.v < 0)
error->all(FLERR,"Variable name for compute global/atom does not exist"); error->all(FLERR,"Variable name {} for compute global/atom index does not exist",
if (input->variable->atomstyle(ivariable) == 0) reference.id);
error->all(FLERR,"Compute global/atom variable is not " if (input->variable->atomstyle(reference.val.v) == 0)
"atom-style variable"); error->all(FLERR,"Compute global/atom index variable {} is not atom-style variable",
reference.id);
} }
for (int i = 0; i < nvalues; i++) { for (auto &val : values) {
if (which[i] == ArgInfo::COMPUTE) { if (val.which == ArgInfo::COMPUTE) {
int icompute = modify->find_compute(ids[i]); val.val.c = modify->get_compute_by_id(val.id);
if (icompute < 0) if (!val.val.c)
error->all(FLERR,"Compute ID for compute global/atom does not exist"); error->all(FLERR,"Compute ID {} for compute global/atom does not exist", val.id);
if (argindex[i] == 0) { if (val.argindex == 0) {
if (!modify->compute[icompute]->vector_flag) if (!val.val.c->vector_flag)
error->all(FLERR,"Compute global/atom compute does not " error->all(FLERR,"Compute ID {} for global/atom compute does not calculate "
"calculate a global vector"); "a global vector", val.id);
} else { } else {
if (!modify->compute[icompute]->array_flag) if (!val.val.c->array_flag)
error->all(FLERR,"Compute global/atom compute does not " error->all(FLERR,"Compute ID {} for global/atom compute does not calculate "
"calculate a global array"); "a global array", val.id);
if (argindex[i] > modify->compute[icompute]->size_array_cols) if (val.argindex > val.val.c->size_array_cols)
error->all(FLERR,"Compute global/atom compute array is " error->all(FLERR,"Compute global/atom compute {} array is accessed out-of-range", val.id);
"accessed out-of-range");
} }
} else if (which[i] == ArgInfo::FIX) { } else if (val.which == ArgInfo::FIX) {
auto ifix = modify->get_fix_by_id(ids[i]); val.val.f = modify->get_fix_by_id(val.id);
if (!ifix) if (!val.val.f) error->all(FLERR,"Fix ID {} for compute global/atom does not exist", val.id);
error->all(FLERR,"Fix ID for compute global/atom does not exist"); if (val.argindex == 0) {
if (argindex[i] == 0) { if (!val.val.f->vector_flag)
if (!ifix->vector_flag) error->all(FLERR,"Fix ID {} for compute global/atom compute does not calculate "
error->all(FLERR,"Compute global/atom fix {} does not calculate a global vector", ids[i]); "a global vector", val.id);
} else { } else {
if (!ifix->array_flag) if (!val.val.f->array_flag)
error->all(FLERR,"Compute global/atom fix {} does not calculate a global array", ids[i]); error->all(FLERR,"Fix ID {} for compute global/atom compute does not calculate "
if (argindex[i] > ifix->size_array_cols) "a global array", val.id);
error->all(FLERR,"Compute global/atom fix {} array is accessed out-of-range", ids[i]); if (val.argindex > val.val.f->size_array_cols)
error->all(FLERR,"Compute global/atom fix {} array is accessed out-of-range", val.id);
} }
} else if (which[i] == ArgInfo::VARIABLE) { } else if (val.which == ArgInfo::VARIABLE) {
int ivariable = input->variable->find(ids[i]); val.val.v = input->variable->find(val.id.c_str());
if (ivariable < 0) if (val.val.v < 0)
error->all(FLERR,"Variable name for compute global/atom does not exist"); error->all(FLERR,"Variable name {} for compute global/atom does not exist", val.id);
if (input->variable->vectorstyle(ivariable) == 0) if (input->variable->vectorstyle(val.val.v) == 0)
error->all(FLERR,"Compute global/atom variable is not vector-style variable"); error->all(FLERR,"Compute global/atom variable {} is not vector-style variable", val.id);
} }
} }
// this compute produces either a peratom vector or array // this compute produces either a peratom vector or array
if (nvalues == 1) size_peratom_cols = 0; if (values.size() == 1) size_peratom_cols = 0;
else size_peratom_cols = nvalues; else size_peratom_cols = values.size();
nmax = maxvector = 0; nmax = maxvector = 0;
vector_atom = nullptr; vector_atom = nullptr;
@ -193,14 +186,6 @@ ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) :
ComputeGlobalAtom::~ComputeGlobalAtom() ComputeGlobalAtom::~ComputeGlobalAtom()
{ {
delete [] idref;
delete [] which;
delete [] argindex;
for (int m = 0; m < nvalues; m++) delete [] ids[m];
delete [] ids;
delete [] value2index;
memory->destroy(indices); memory->destroy(indices);
memory->destroy(varatom); memory->destroy(varatom);
memory->destroy(vecglobal); memory->destroy(vecglobal);
@ -212,44 +197,38 @@ ComputeGlobalAtom::~ComputeGlobalAtom()
void ComputeGlobalAtom::init() void ComputeGlobalAtom::init()
{ {
// set indices of all computes,fixes,variables // store references of all computes, fixes, or variables
if (whichref == ArgInfo::COMPUTE) { if (reference.which == ArgInfo::COMPUTE) {
int icompute = modify->find_compute(idref); reference.val.c = modify->get_compute_by_id(reference.id);
if (icompute < 0) if (!reference.val.c)
error->all(FLERR,"Compute ID for compute global/atom does not exist"); error->all(FLERR,"Compute ID {} for compute global/atom index does not exist", reference.id);
ref2index = icompute; } else if (reference.which == ArgInfo::FIX) {
} else if (whichref == ArgInfo::FIX) { reference.val.f = modify->get_fix_by_id(reference.id);
int ifix = modify->find_fix(idref); if (reference.val.f)
if (ifix < 0) error->all(FLERR,"Fix ID {} for compute global/atom index does not exist", reference.id);
error->all(FLERR,"Fix ID for compute global/atom does not exist"); } else if (reference.which == ArgInfo::VARIABLE) {
ref2index = ifix; reference.val.v = input->variable->find(reference.id.c_str());
} else if (whichref == ArgInfo::VARIABLE) { if (reference.val.v < 0)
int ivariable = input->variable->find(idref); error->all(FLERR,"Variable name {} for compute global/atom index does not exist",
if (ivariable < 0) reference.id);
error->all(FLERR,"Variable name for compute global/atom does not exist");
ref2index = ivariable;
} }
for (int m = 0; m < nvalues; m++) { for (auto &val : values) {
if (which[m] == ArgInfo::COMPUTE) { if (val.which == ArgInfo::COMPUTE) {
int icompute = modify->find_compute(ids[m]); val.val.c = modify->get_compute_by_id(val.id);
if (icompute < 0) if (!val.val.c)
error->all(FLERR,"Compute ID for compute global/atom does not exist"); error->all(FLERR,"Compute ID {} for compute global/atom does not exist", val.id);
value2index[m] = icompute;
} else if (which[m] == ArgInfo::FIX) { } else if (val.which == ArgInfo::FIX) {
int ifix = modify->find_fix(ids[m]); val.val.f = modify->get_fix_by_id(val.id);
if (ifix < 0) if (!val.val.f)
error->all(FLERR,"Fix ID for compute global/atom does not exist"); error->all(FLERR,"Fix ID {} for compute global/atom does not exist", val.id);
value2index[m] = ifix;
} else if (which[m] == ArgInfo::VARIABLE) { } else if (val.which == ArgInfo::VARIABLE) {
int ivariable = input->variable->find(ids[m]); val.val.v = input->variable->find(val.id.c_str());
if (ivariable < 0) if (val.val.v < 0)
error->all(FLERR,"Variable name for compute global/atom " error->all(FLERR,"Variable name {} for compute global/atom does not exist", val.id);
"does not exist");
value2index[m] = ivariable;
} }
} }
} }
@ -258,8 +237,6 @@ void ComputeGlobalAtom::init()
void ComputeGlobalAtom::compute_peratom() void ComputeGlobalAtom::compute_peratom()
{ {
int i,j;
invoked_peratom = update->ntimestep; invoked_peratom = update->ntimestep;
// grow indices and output vector or array if necessary // grow indices and output vector or array if necessary
@ -268,16 +245,16 @@ void ComputeGlobalAtom::compute_peratom()
nmax = atom->nmax; nmax = atom->nmax;
memory->destroy(indices); memory->destroy(indices);
memory->create(indices,nmax,"global/atom:indices"); memory->create(indices,nmax,"global/atom:indices");
if (whichref == ArgInfo::VARIABLE) { if (reference.which == ArgInfo::VARIABLE) {
memory->destroy(varatom); memory->destroy(varatom);
memory->create(varatom,nmax,"global/atom:varatom"); memory->create(varatom,nmax,"global/atom:varatom");
} }
if (nvalues == 1) { if (values.size() == 1) {
memory->destroy(vector_atom); memory->destroy(vector_atom);
memory->create(vector_atom,nmax,"global/atom:vector_atom"); memory->create(vector_atom,nmax,"global/atom:vector_atom");
} else { } else {
memory->destroy(array_atom); memory->destroy(array_atom);
memory->create(array_atom,nmax,nvalues,"global/atom:array_atom"); memory->create(array_atom,nmax,values.size(),"global/atom:array_atom");
} }
} }
@ -286,81 +263,79 @@ void ComputeGlobalAtom::compute_peratom()
// indices are decremented from 1 to N -> 0 to N-1 // indices are decremented from 1 to N -> 0 to N-1
int *mask = atom->mask; int *mask = atom->mask;
int nlocal = atom->nlocal; const int nlocal = atom->nlocal;
if (whichref == ArgInfo::COMPUTE) { if (reference.which == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[ref2index];
if (!(compute->invoked_flag & Compute::INVOKED_PERATOM)) { if (!(reference.val.c->invoked_flag & Compute::INVOKED_PERATOM)) {
compute->compute_peratom(); reference.val.c->compute_peratom();
compute->invoked_flag |= Compute::INVOKED_PERATOM; reference.val.c->invoked_flag |= Compute::INVOKED_PERATOM;
} }
if (indexref == 0) { if (reference.argindex == 0) {
double *compute_vector = compute->vector_atom; double *compute_vector = reference.val.c->vector_atom;
for (i = 0; i < nlocal; i++) for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) if (mask[i] & groupbit)
indices[i] = static_cast<int> (compute_vector[i]) - 1; indices[i] = static_cast<int> (compute_vector[i]) - 1;
} else { } else {
double **compute_array = compute->array_atom; double **compute_array = reference.val.c->array_atom;
int im1 = indexref - 1; int im1 = reference.argindex - 1;
for (i = 0; i < nlocal; i++) for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) if (mask[i] & groupbit)
indices[i] = static_cast<int> (compute_array[i][im1]) - 1; indices[i] = static_cast<int> (compute_array[i][im1]) - 1;
} }
} else if (whichref == ArgInfo::FIX) { } else if (reference.which == ArgInfo::FIX) {
if (update->ntimestep % modify->fix[ref2index]->peratom_freq) if (update->ntimestep % reference.val.f->peratom_freq)
error->all(FLERR,"Fix used in compute global/atom not " error->all(FLERR,"Fix {} used in compute global/atom not computed at compatible time",
"computed at compatible time"); reference.id);
Fix *fix = modify->fix[ref2index];
if (indexref == 0) { if (reference.argindex == 0) {
double *fix_vector = fix->vector_atom; double *fix_vector = reference.val.f->vector_atom;
for (i = 0; i < nlocal; i++) for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) if (mask[i] & groupbit)
indices[i] = static_cast<int> (fix_vector[i]) - 1; indices[i] = static_cast<int> (fix_vector[i]) - 1;
} else { } else {
double **fix_array = fix->array_atom; double **fix_array = reference.val.f->array_atom;
int im1 = indexref - 1; int im1 = reference.argindex - 1;
for (i = 0; i < nlocal; i++) for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) if (mask[i] & groupbit)
indices[i] = static_cast<int> (fix_array[i][im1]) - 1; indices[i] = static_cast<int> (fix_array[i][im1]) - 1;
} }
} else if (whichref == ArgInfo::VARIABLE) { } else if (reference.which == ArgInfo::VARIABLE) {
input->variable->compute_atom(ref2index,igroup,varatom,1,0); input->variable->compute_atom(reference.val.v, igroup, varatom, 1, 0);
for (i = 0; i < nlocal; i++) for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) if (mask[i] & groupbit)
indices[i] = static_cast<int> (varatom[i]) - 1; indices[i] = static_cast<int> (varatom[i]) - 1;
} }
// loop over values to fill output vector or array // loop over values to fill output vector or array
for (int m = 0; m < nvalues; m++) { int m = 0;
for (auto &val : values) {
// output = vector // output = vector
if (argindex[m] == 0) { if (val.argindex == 0) {
int vmax; int vmax;
double *source; double *source;
if (which[m] == ArgInfo::COMPUTE) { if (val.which == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[value2index[m]];
if (!(compute->invoked_flag & Compute::INVOKED_VECTOR)) { if (!(val.val.c->invoked_flag & Compute::INVOKED_VECTOR)) {
compute->compute_vector(); val.val.c->compute_vector();
compute->invoked_flag |= Compute::INVOKED_VECTOR; val.val.c->invoked_flag |= Compute::INVOKED_VECTOR;
} }
source = compute->vector; source = val.val.c->vector;
vmax = compute->size_vector; vmax = val.val.c->size_vector;
} else if (which[m] == ArgInfo::FIX) { } else if (val.which == ArgInfo::FIX) {
if (update->ntimestep % modify->fix[value2index[m]]->peratom_freq) if (update->ntimestep % val.val.f->peratom_freq)
error->all(FLERR,"Fix used in compute global/atom not computed at compatible time"); error->all(FLERR,"Fix {} used in compute global/atom not computed at compatible time",
Fix *fix = modify->fix[value2index[m]]; val.id);
vmax = fix->size_vector; vmax = reference.val.f->size_vector;
if (vmax > maxvector) { if (vmax > maxvector) {
maxvector = vmax; maxvector = vmax;
@ -368,28 +343,28 @@ void ComputeGlobalAtom::compute_peratom()
memory->create(vecglobal,maxvector,"global/atom:vecglobal"); memory->create(vecglobal,maxvector,"global/atom:vecglobal");
} }
for (i = 0; i < vmax; i++) for (int i = 0; i < vmax; i++)
vecglobal[i] = fix->compute_vector(i); vecglobal[i] = val.val.f->compute_vector(i);
source = vecglobal; source = vecglobal;
} else if (which[m] == ArgInfo::VARIABLE) { } else if (val.which == ArgInfo::VARIABLE) {
vmax = input->variable->compute_vector(value2index[m],&source); vmax = input->variable->compute_vector(val.val.v, &source);
} }
if (nvalues == 1) { if (values.size() == 1) {
for (i = 0; i < nlocal; i++) { for (int i = 0; i < nlocal; i++) {
vector_atom[i] = 0.0; vector_atom[i] = 0.0;
if (mask[i] & groupbit) { if (mask[i] & groupbit) {
j = indices[i]; int j = indices[i];
if (j >= 0 && j < vmax) vector_atom[i] = source[j]; if (j >= 0 && j < vmax) vector_atom[i] = source[j];
} }
} }
} else { } else {
for (i = 0; i < nlocal; i++) { for (int i = 0; i < nlocal; i++) {
array_atom[i][m] = 0.0; array_atom[i][m] = 0.0;
if (mask[i] & groupbit) { if (mask[i] & groupbit) {
j = indices[i]; int j = indices[i];
if (j >= 0 && j < vmax) array_atom[i][m] = source[j]; if (j >= 0 && j < vmax) array_atom[i][m] = source[j];
} }
} }
@ -400,18 +375,17 @@ void ComputeGlobalAtom::compute_peratom()
} else { } else {
int vmax; int vmax;
double *source; double *source;
int col = argindex[m] - 1; int col = val.argindex - 1;
if (which[m] == ArgInfo::COMPUTE) { if (val.which == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[value2index[m]];
if (!(compute->invoked_flag & Compute::INVOKED_ARRAY)) { if (!(val.val.c->invoked_flag & Compute::INVOKED_ARRAY)) {
compute->compute_array(); val.val.c->compute_array();
compute->invoked_flag |= Compute::INVOKED_ARRAY; val.val.c->invoked_flag |= Compute::INVOKED_ARRAY;
} }
double **compute_array = compute->array; double **compute_array = val.val.c->array;
vmax = compute->size_array_rows; vmax = val.val.c->size_array_rows;
if (vmax > maxvector) { if (vmax > maxvector) {
maxvector = vmax; maxvector = vmax;
@ -419,17 +393,16 @@ void ComputeGlobalAtom::compute_peratom()
memory->create(vecglobal,maxvector,"global/atom:vecglobal"); memory->create(vecglobal,maxvector,"global/atom:vecglobal");
} }
for (i = 0; i < vmax; i++) for (int i = 0; i < vmax; i++)
vecglobal[i] = compute_array[i][col]; vecglobal[i] = compute_array[i][col];
source = vecglobal; source = vecglobal;
} else if (which[m] == ArgInfo::FIX) { } else if (val.which == ArgInfo::FIX) {
if (update->ntimestep % modify->fix[value2index[m]]->peratom_freq) if (update->ntimestep % val.val.f->peratom_freq)
error->all(FLERR,"Fix used in compute global/atom not " error->all(FLERR,"Fix {} used in compute global/atom not computed at compatible time",
"computed at compatible time"); val.id);
Fix *fix = modify->fix[value2index[m]]; vmax = val.val.f->size_array_rows;
vmax = fix->size_array_rows;
if (vmax > maxvector) { if (vmax > maxvector) {
maxvector = vmax; maxvector = vmax;
@ -437,34 +410,35 @@ void ComputeGlobalAtom::compute_peratom()
memory->create(vecglobal,maxvector,"global/atom:vecglobal"); memory->create(vecglobal,maxvector,"global/atom:vecglobal");
} }
for (i = 0; i < vmax; i++) for (int i = 0; i < vmax; i++)
vecglobal[i] = fix->compute_array(i,col); vecglobal[i] = val.val.f->compute_array(i,col);
source = vecglobal; source = vecglobal;
} else if (which[m] == ArgInfo::VARIABLE) { } else if (val.which == ArgInfo::VARIABLE) {
vmax = input->variable->compute_vector(value2index[m],&source); vmax = input->variable->compute_vector(val.val.v, &source);
} }
if (nvalues == 1) { if (values.size() == 1) {
for (i = 0; i < nlocal; i++) { for (int i = 0; i < nlocal; i++) {
vector_atom[i] = 0.0; vector_atom[i] = 0.0;
if (mask[i] & groupbit) { if (mask[i] & groupbit) {
j = indices[i]; int j = indices[i];
if (j >= 0 && j < vmax) vector_atom[i] = source[j]; if (j >= 0 && j < vmax) vector_atom[i] = source[j];
} }
} }
} else { } else {
for (i = 0; i < nlocal; i++) { for (int i = 0; i < nlocal; i++) {
array_atom[i][m] = 0.0; array_atom[i][m] = 0.0;
if (mask[i] & groupbit) { if (mask[i] & groupbit) {
j = indices[i]; int j = indices[i];
if (j >= 0 && j < vmax) array_atom[i][m] = source[j]; if (j >= 0 && j < vmax) array_atom[i][m] = source[j];
} }
} }
} }
} }
} }
++m;
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -473,7 +447,7 @@ void ComputeGlobalAtom::compute_peratom()
double ComputeGlobalAtom::memory_usage() double ComputeGlobalAtom::memory_usage()
{ {
double bytes = (double)nmax*nvalues * sizeof(double); double bytes = (double)nmax*values.size() * sizeof(double);
bytes += (double)nmax * sizeof(int); // indices bytes += (double)nmax * sizeof(int); // indices
if (varatom) bytes += (double)nmax * sizeof(double); // varatom if (varatom) bytes += (double)nmax * sizeof(double); // varatom
bytes += (double)maxvector * sizeof(double); // vecglobal bytes += (double)maxvector * sizeof(double); // vecglobal

View File

@ -33,12 +33,18 @@ class ComputeGlobalAtom : public Compute {
double memory_usage() override; double memory_usage() override;
protected: protected:
int whichref, indexref, ref2index; struct value_t {
char *idref; int which;
int argindex;
int nvalues; std::string id;
int *which, *argindex, *value2index; union {
char **ids; class Compute *c;
class Fix *f;
int v;
} val;
};
std::vector<value_t> values;
value_t reference;
int nmax, maxvector; int nmax, maxvector;
int *indices; int *indices;

View File

@ -648,7 +648,6 @@ TEST_F(AtomStyleTest, atomic)
END_HIDE_OUTPUT(); END_HIDE_OUTPUT();
ASSERT_EQ(lmp->atom->map_tag_max, 16); ASSERT_EQ(lmp->atom->map_tag_max, 16);
x = lmp->atom->x; x = lmp->atom->x;
tag = lmp->atom->tag;
EXPECT_NEAR(x[GETIDX(1)][0], -2.0, EPSILON); EXPECT_NEAR(x[GETIDX(1)][0], -2.0, EPSILON);
EXPECT_NEAR(x[GETIDX(1)][1], 2.0, EPSILON); EXPECT_NEAR(x[GETIDX(1)][1], 2.0, EPSILON);
EXPECT_NEAR(x[GETIDX(1)][2], 0.1, EPSILON); EXPECT_NEAR(x[GETIDX(1)][2], 0.1, EPSILON);