ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

refactor fix ave/time to use array of struct

Browse Source
This commit is contained in:
Axel Kohlmeyer
2022-10-02 14:48:53 -04:00
parent 39ff5c5722
commit d7ff6337ed
2 changed files with 313 additions and 345 deletions
Download Patch File Download Diff File Expand all files Collapse all files

633
src/fix_ave_time.cpp
View File

@ -19,6 +19,7 @@
#include "fix_ave_time.h"
#include "arg_info.h"
#include "comm.h"
#include "compute.h"
#include "error.h"
#include "input.h"
@ -32,22 +33,18 @@
using namespace LAMMPS_NS;
using namespace FixConst;
enum{ONE,RUNNING,WINDOW};
enum{SCALAR,VECTOR};
enum{ ONE, RUNNING, WINDOW };
enum{ SCALAR, VECTOR };
/* ---------------------------------------------------------------------- */
FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
Fix(lmp, narg, arg),
nvalues(0), which(nullptr), argindex(nullptr), value2index(nullptr),
offcol(nullptr), varlen(nullptr), ids(nullptr),
fp(nullptr), offlist(nullptr), format(nullptr), format_user(nullptr),
nvalues(0), fp(nullptr), offlist(nullptr), format(nullptr), format_user(nullptr),
vector(nullptr), vector_total(nullptr), vector_list(nullptr),
column(nullptr), array(nullptr), array_total(nullptr), array_list(nullptr)
{
if (narg < 7) error->all(FLERR,"Illegal fix ave/time command");
MPI_Comm_rank(world,&me);
if (narg < 7) utils::missing_cmd_args(FLERR, "fix ave/time", error);
nevery = utils::inumeric(FLERR,arg[3],false,lmp);
nrepeat = utils::inumeric(FLERR,arg[4],false,lmp);
@ -62,18 +59,17 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
// then read options so know mode = SCALAR/VECTOR before re-reading values
nvalues = 0;
int iarg = 6;
while (iarg < narg) {
if ((strncmp(arg[iarg],"c_",2) == 0) ||
(strncmp(arg[iarg],"f_",2) == 0) ||
(strncmp(arg[iarg],"v_",2) == 0)) {
if (utils::strmatch(arg[iarg],"^[cfv]_")) {
nvalues++;
iarg++;
} else break;
}
if (nvalues == 0)
error->all(FLERR,"No values from computes, fixes, or variables used in fix ave/time command");
if (nvalues == 0) error->all(FLERR,"No values in fix ave/time command");
// parse optional keywords
options(iarg,narg,arg);
@ -83,7 +79,6 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
int expand = 0;
char **earg;
nvalues = utils::expand_args(FLERR,nvalues,&arg[6],mode,earg,lmp);
keyword.resize(nvalues);
key2col.clear();
if (earg != &arg[6]) expand = 1;
@ -91,132 +86,126 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
// parse values
which = new int[nvalues];
argindex = new int[nvalues];
value2index = new int[nvalues];
offcol = new int[nvalues];
varlen = new int[nvalues];
ids = new char*[nvalues];
values.clear();
for (int i = 0; i < nvalues; i++) {
ArgInfo argi(arg[i]);
keyword[i] = arg[i];
value_t val;
val.keyword = arg[i];
key2col[arg[i]] = i;
if ((argi.get_type() == ArgInfo::NONE)
|| (argi.get_type() == ArgInfo::UNKNOWN)
|| (argi.get_dim() > 1))
error->all(FLERR,"Invalid fix ave/time command");
error->all(FLERR,"Invalid fix ave/time argument: {}", arg[i]);
which[i] = argi.get_type();
argindex[i] = argi.get_index1();
ids[i] = argi.copy_name();
val.which = argi.get_type();
val.argindex = argi.get_index1();
val.varlen = 0;
val.offcol = 0;
val.id = argi.get_name();
val.val.c = nullptr;
values.push_back(val);
}
if (nvalues != (int)values.size())
error->all(FLERR, "Could not parse value data consistently for fix ave/time");
// set off columns now that nvalues is finalized
for (int i = 0; i < nvalues; i++) offcol[i] = 0;
for (int i = 0; i < noff; i++) {
if (offlist[i] < 1 || offlist[i] > nvalues)
error->all(FLERR,"Invalid fix ave/time off column");
offcol[offlist[i]-1] = 1;
error->all(FLERR,"Invalid fix ave/time off column: {}", offlist[i]);
values[offlist[i]-1].offcol = 1;
}
// setup and error check
// for fix inputs, check that fix frequency is acceptable
// set variable_length if any compute is variable length
if (nevery <= 0 || nrepeat <= 0 || nfreq <= 0)
error->all(FLERR,"Illegal fix ave/time command");
if (nevery <= 0)
error->all(FLERR,"Illegal fix ave/time nevery value: {}", nevery);
if (nrepeat <= 0)
error->all(FLERR,"Illegal fix ave/time nrepeat value: {}", nrepeat);
if (nfreq <= 0)
error->all(FLERR,"Illegal fix ave/time nfreq value: {}", nfreq);
if (nfreq % nevery || nrepeat*nevery > nfreq)
error->all(FLERR,"Illegal fix ave/time command");
error->all(FLERR,"Inconsistent fix ave/time nevery/nrepeat/nfreq values");
if (ave != RUNNING && overwrite)
error->all(FLERR,"Illegal fix ave/time command");
error->all(FLERR,"Fix ave/time overwrite keyword requires ave running setting");
for (int i = 0; i < nvalues; i++) {
varlen[i] = 0;
for (auto &val : values) {
if (which[i] == ArgInfo::COMPUTE && mode == SCALAR) {
int icompute = modify->find_compute(ids[i]);
if (icompute < 0)
error->all(FLERR,"Compute ID for fix ave/time does not exist");
if (argindex[i] == 0 && modify->compute[icompute]->scalar_flag == 0)
error->all(FLERR,"Fix ave/time compute does not calculate a scalar");
if (argindex[i] && modify->compute[icompute]->vector_flag == 0)
error->all(FLERR,"Fix ave/time compute does not calculate a vector");
if (argindex[i] && argindex[i] > modify->compute[icompute]->size_vector &&
modify->compute[icompute]->size_vector_variable == 0)
error->all(FLERR,
"Fix ave/time compute vector is accessed out-of-range");
if (argindex[i] && modify->compute[icompute]->size_vector_variable)
varlen[i] = 1;
if ((val.which == ArgInfo::COMPUTE) && (mode == SCALAR)) {
val.val.c = modify->get_compute_by_id(val.id);
if (!val.val.c) error->all(FLERR,"Compute ID {} for fix ave/time does not exist", val.id);
if (val.argindex == 0 && (val.val.c->scalar_flag == 0))
error->all(FLERR,"Fix ave/time compute {} does not calculate a scalar", val.id);
if (val.argindex && (val.val.c->vector_flag == 0))
error->all(FLERR,"Fix ave/time compute {} does not calculate a vector", val.id);
if (val.argindex && (val.argindex > val.val.c->size_vector) &&
(val.val.c->size_vector_variable == 0))
error->all(FLERR, "Fix ave/time compute {} vector is accessed out-of-range", val.id);
if (val.argindex && val.val.c->size_vector_variable) val.varlen = 1;
} else if (which[i] == ArgInfo::COMPUTE && mode == VECTOR) {
int icompute = modify->find_compute(ids[i]);
if (icompute < 0)
error->all(FLERR,"Compute ID for fix ave/time does not exist");
if (argindex[i] == 0 && modify->compute[icompute]->vector_flag == 0)
error->all(FLERR,"Fix ave/time compute does not calculate a vector");
if (argindex[i] && modify->compute[icompute]->array_flag == 0)
error->all(FLERR,"Fix ave/time compute does not calculate an array");
if (argindex[i] &&
argindex[i] > modify->compute[icompute]->size_array_cols)
error->all(FLERR,"Fix ave/time compute array is accessed out-of-range");
if (argindex[i] == 0 && modify->compute[icompute]->size_vector_variable)
varlen[i] = 1;
if (argindex[i] && modify->compute[icompute]->size_array_rows_variable)
varlen[i] = 1;
} else if ((val.which == ArgInfo::COMPUTE) && (mode == VECTOR)) {
val.val.c = modify->get_compute_by_id(val.id);
if (!val.val.c) error->all(FLERR,"Compute ID {} for fix ave/time does not exist", val.id);
if ((val.argindex == 0) && (val.val.c->vector_flag == 0))
error->all(FLERR,"Fix ave/time compute {} does not calculate a vector", val.id);
if (val.argindex && (val.val.c->array_flag == 0))
error->all(FLERR,"Fix ave/time compute {} does not calculate an array", val.id);
if (val.argindex && (val.argindex > val.val.c->size_array_cols))
error->all(FLERR,"Fix ave/time compute {} array is accessed out-of-range", val.id);
if ((val.argindex == 0) && (val.val.c->size_vector_variable)) val.varlen = 1;
if (val.argindex && (val.val.c->size_array_rows_variable)) val.varlen = 1;
} else if (which[i] == ArgInfo::FIX && mode == SCALAR) {
int ifix = modify->find_fix(ids[i]);
if (ifix < 0)
error->all(FLERR,"Fix ID for fix ave/time does not exist");
if (argindex[i] == 0 && modify->fix[ifix]->scalar_flag == 0)
error->all(FLERR,"Fix ave/time fix does not calculate a scalar");
if (argindex[i] && modify->fix[ifix]->vector_flag == 0)
error->all(FLERR,"Fix ave/time fix does not calculate a vector");
if (argindex[i] && modify->fix[ifix]->size_vector_variable)
error->all(FLERR,"Fix ave/time fix vector cannot be variable length");
if (argindex[i] && argindex[i] > modify->fix[ifix]->size_vector)
error->all(FLERR,"Fix ave/time fix vector is accessed out-of-range");
if (nevery % modify->fix[ifix]->global_freq)
error->all(FLERR,
"Fix for fix ave/time not computed at compatible time");
} else if ((val.which == ArgInfo::FIX) && (mode == SCALAR)) {
val.val.f = modify->get_fix_by_id(val.id);
if (!val.val.f) error->all(FLERR,"Fix ID {} for fix ave/time does not exist", val.id);
if ((val.argindex == 0) && (val.val.f->scalar_flag == 0))
error->all(FLERR,"Fix ave/time fix {} does not calculate a scalar", val.id);
if (val.argindex && (val.val.f->vector_flag == 0))
error->all(FLERR,"Fix ave/time fix {} does not calculate a vector", val.id);
if (val.argindex && (val.val.f->size_vector_variable))
error->all(FLERR,"Fix ave/time fix {} vector cannot be variable length", val.id);
if (val.argindex && (val.argindex > val.val.f->size_vector))
error->all(FLERR,"Fix ave/time fix {} vector is accessed out-of-range", val.id);
if (nevery % val.val.f->global_freq)
error->all(FLERR, "Fix {} for fix ave/time not computed at compatible time", val.id);
} else if (which[i] == ArgInfo::FIX && mode == VECTOR) {
int ifix = modify->find_fix(ids[i]);
if (ifix < 0)
error->all(FLERR,"Fix ID for fix ave/time does not exist");
if (argindex[i] == 0 && modify->fix[ifix]->vector_flag == 0)
error->all(FLERR,"Fix ave/time fix does not calculate a vector");
if (argindex[i] && modify->fix[ifix]->array_flag == 0)
error->all(FLERR,"Fix ave/time fix does not calculate an array");
if (argindex[i] && modify->fix[ifix]->size_array_rows_variable)
error->all(FLERR,"Fix ave/time fix array cannot be variable length");
if (argindex[i] && argindex[i] > modify->fix[ifix]->size_array_cols)
error->all(FLERR,"Fix ave/time fix array is accessed out-of-range");
if (nevery % modify->fix[ifix]->global_freq)
error->all(FLERR,
"Fix for fix ave/time not computed at compatible time");
} else if ((val.which == ArgInfo::FIX) && (mode == VECTOR)) {
val.val.f = modify->get_fix_by_id(val.id);
if (!val.val.f) error->all(FLERR,"Fix ID {} for fix ave/time does not exist", val.id);
if ((val.argindex == 0) && (val.val.f->vector_flag == 0))
error->all(FLERR,"Fix ave/time fix {} does not calculate a vector", val.id);
if (val.argindex && (val.val.f->array_flag == 0))
error->all(FLERR,"Fix ave/time fix {} does not calculate an array", val.id);
if (val.argindex && (val.val.f->size_array_rows_variable))
error->all(FLERR,"Fix ave/time fix {} array cannot be variable length", val.id);
if (val.argindex && (val.argindex > val.val.f->size_array_cols))
error->all(FLERR,"Fix ave/time fix {} array is accessed out-of-range", val.id);
if (nevery % val.val.f->global_freq)
error->all(FLERR, "Fix {} for fix ave/time not computed at compatible time", val.id);
} else if (which[i] == ArgInfo::VARIABLE && mode == SCALAR) {
int ivariable = input->variable->find(ids[i]);
} else if ((val.which == ArgInfo::VARIABLE) && (mode == SCALAR)) {
int ivariable = input->variable->find(val.id.c_str());
if (ivariable < 0)
error->all(FLERR,"Variable name for fix ave/time does not exist");
if (argindex[i] == 0 && input->variable->equalstyle(ivariable) == 0)
error->all(FLERR,"Fix ave/time variable is not equal-style variable");
if (argindex[i] && input->variable->vectorstyle(ivariable) == 0)
error->all(FLERR,"Fix ave/time variable is not vector-style variable");
error->all(FLERR,"Variable name {} for fix ave/time does not exist", val.id);
if ((val.argindex == 0) && (input->variable->equalstyle(ivariable) == 0))
error->all(FLERR,"Fix ave/time variable {} is not equal-style variable", val.id);
if ((val.argindex) && (input->variable->vectorstyle(ivariable) == 0))
error->all(FLERR,"Fix ave/time variable {} is not vector-style variable", val.id);
} else if (which[i] == ArgInfo::VARIABLE && mode == VECTOR) {
int ivariable = input->variable->find(ids[i]);
} else if ((val.which == ArgInfo::VARIABLE) && (mode == VECTOR)) {
int ivariable = input->variable->find(val.id.c_str());
if (ivariable < 0)
error->all(FLERR,"Variable name for fix ave/time does not exist");
if (argindex[i] == 0 && input->variable->vectorstyle(ivariable) == 0)
error->all(FLERR,"Fix ave/time variable is not vector-style variable");
if (argindex[i])
error->all(FLERR,"Fix ave/time mode vector variable cannot be indexed");
varlen[i] = 1;
error->all(FLERR,"Variable name {} for fix ave/time does not exist", val.id);
if ((val.argindex == 0) && (input->variable->vectorstyle(ivariable) == 0))
error->all(FLERR,"Fix ave/time variable {} is not vector-style variable", val.id);
if (val.argindex)
error->all(FLERR,"Fix ave/time mode vector variable {} cannot be indexed", val.id);
val.varlen = 1;
}
}
@ -225,9 +214,9 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
all_variable_length = 1;
any_variable_length = 0;
for (int i = 0; i < nvalues; i++) {
if (varlen[i] == 0) all_variable_length = 0;
if (varlen[i]) any_variable_length = 1;
for (auto &val : values) {
if (val.varlen == 0) all_variable_length = 0;
if (val.varlen) any_variable_length = 1;
}
// if VECTOR mode, check that all columns are same length
@ -245,19 +234,18 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
// only if nrepeat > 1 or ave = RUNNING/WINDOW,
// so that locking spans multiple timesteps
if (any_variable_length &&
(nrepeat > 1 || ave == RUNNING || ave == WINDOW)) {
for (int i = 0; i < nvalues; i++)
if (varlen[i] && which[i] == ArgInfo::COMPUTE)
modify->get_compute_by_id(ids[i])->lock_enable();
lockforever = 0;
if (any_variable_length && ((nrepeat > 1) || (ave == RUNNING) || (ave == WINDOW))) {
for (auto &val : values) {
if (val.varlen && val.which == ArgInfo::COMPUTE) val.val.c->lock_enable();
lockforever = 0;
}
}
// print file comment lines
// for mode = VECTOR, cannot use arg to print
// since array args may have been expanded to multiple vectors
if (fp && me == 0) {
if (fp && comm->me == 0) {
clearerr(fp);
if (title1) fprintf(fp,"%s\n",title1);
else fprintf(fp,"# Time-averaged data for fix %s\n",id);
@ -300,8 +288,7 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
if (mode == SCALAR) {
vector = new double[nvalues];
vector_total = new double[nvalues];
if (ave == WINDOW)
memory->create(vector_list,nwindow,nvalues,"ave/time:vector_list");
if (ave == WINDOW) memory->create(vector_list,nwindow,nvalues,"ave/time:vector_list");
} else allocate_arrays();
// this fix produces either a global scalar or vector or array
@ -314,17 +301,16 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
if (mode == SCALAR) {
if (nvalues == 1) {
scalar_flag = 1;
if (which[0] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[modify->find_compute(ids[0])];
if (argindex[0] == 0) extscalar = compute->extscalar;
else if (compute->extvector >= 0) extscalar = compute->extvector;
else extscalar = compute->extlist[argindex[0]-1];
} else if (which[0] == ArgInfo::FIX) {
Fix *fix = modify->fix[modify->find_fix(ids[0])];
if (argindex[0] == 0) extscalar = fix->extscalar;
else if (fix->extvector >= 0) extscalar = fix->extvector;
else extscalar = fix->extlist[argindex[0]-1];
} else if (which[0] == ArgInfo::VARIABLE) {
auto &val = values[0];
if (val.which == ArgInfo::COMPUTE) {
if (val.argindex == 0) extscalar = val.val.c->extscalar;
else if (val.val.c->extvector >= 0) extscalar = val.val.c->extvector;
else extscalar = val.val.c->extlist[val.argindex-1];
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0) extscalar = val.val.f->extscalar;
else if (val.val.f->extvector >= 0) extscalar = val.val.f->extvector;
else extscalar = val.val.f->extlist[val.argindex-1];
} else if (val.which == ArgInfo::VARIABLE) {
extscalar = 0;
}
@ -333,47 +319,46 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
size_vector = nrows = nvalues;
extvector = -1;
extlist = new int[nvalues];
for (int i = 0; i < nvalues; i++) {
if (which[i] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[modify->find_compute(ids[i])];
if (argindex[i] == 0) extlist[i] = compute->extscalar;
else if (compute->extvector >= 0) extlist[i] = compute->extvector;
else extlist[i] = compute->extlist[argindex[i]-1];
} else if (which[i] == ArgInfo::FIX) {
Fix *fix = modify->fix[modify->find_fix(ids[i])];
if (argindex[i] == 0) extlist[i] = fix->extscalar;
else if (fix->extvector >= 0) extlist[i] = fix->extvector;
else extlist[i] = fix->extlist[argindex[i]-1];
} else if (which[i] == ArgInfo::VARIABLE) {
int i = 0;
for (auto &val : values) {
if (val.which == ArgInfo::COMPUTE) {
if (val.argindex == 0) extlist[i] = val.val.c->extscalar;
else if (val.val.c->extvector >= 0) extlist[i] = val.val.c->extvector;
else extlist[i] = val.val.c->extlist[val.argindex-1];
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0) extlist[i] = val.val.f->extscalar;
else if (val.val.f->extvector >= 0) extlist[i] = val.val.f->extvector;
else extlist[i] = val.val.f->extlist[val.argindex-1];
} else if (val.which == ArgInfo::VARIABLE) {
extlist[i] = 0;
}
++i;
}
}
} else {
if (nvalues == 1) {
auto &val = values[0];
vector_flag = 1;
size_vector = nrows;
if (all_variable_length) size_vector_variable = 1;
if (which[0] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[modify->find_compute(ids[0])];
if (argindex[0] == 0) {
extvector = compute->extvector;
if (val.which == ArgInfo::COMPUTE) {
if (val.argindex == 0) {
extvector = val.val.c->extvector;
if (extvector == -1) {
extlist = new int[nrows];
for (int i = 0; i < nrows; i++) extlist[i] = compute->extlist[i];
for (int i = 0; i < nrows; i++) extlist[i] = val.val.c->extlist[i];
}
} else extvector = compute->extarray;
} else if (which[0] == ArgInfo::FIX) {
Fix *fix = modify->fix[modify->find_fix(ids[0])];
if (argindex[0] == 0) {
extvector = fix->extvector;
} else extvector = val.val.c->extarray;
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0) {
extvector = val.val.f->extvector;
if (extvector == -1) {
extlist = new int[nrows];
for (int i = 0; i < nrows; i++) extlist[i] = fix->extlist[i];
for (int i = 0; i < nrows; i++) extlist[i] = val.val.f->extlist[i];
}
} else extvector = fix->extarray;
} else if (which[0] == ArgInfo::VARIABLE) {
} else extvector = val.val.f->extarray;
} else if (val.which == ArgInfo::VARIABLE) {
extlist = new int[nrows];
for (int i = 0; i < nrows; i++) extlist[i] = 0;
}
@ -383,26 +368,25 @@ FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
size_array_rows = nrows;
size_array_cols = nvalues;
if (all_variable_length) size_array_rows_variable = 1;
int value;
for (int i = 0; i < nvalues; i++) {
if (which[i] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[modify->find_compute(ids[i])];
if (argindex[i] == 0) value = compute->extvector;
else value = compute->extarray;
} else if (which[i] == ArgInfo::FIX) {
Fix *fix = modify->fix[modify->find_fix(ids[i])];
if (argindex[i] == 0) value = fix->extvector;
else value = fix->extarray;
} else if (which[i] == ArgInfo::VARIABLE) {
value = 0;
int extvalue = 0;
extarray = -2;
for (auto &val : values) {
if (val.which == ArgInfo::COMPUTE) {
if (val.argindex == 0) extvalue = val.val.c->extvector;
else extvalue = val.val.c->extarray;
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0) extvalue = val.val.f->extvector;
else extvalue = val.val.f->extarray;
} else if (val.which == ArgInfo::VARIABLE) {
extvalue = 0;
}
if (value == -1)
error->all(FLERR,"Fix ave/time cannot set output array "
"intensive/extensive from these inputs");
if (i == 0) extarray = value;
else if (value != extarray)
error->all(FLERR,"Fix ave/time cannot set output array "
"intensive/extensive from these inputs");
if (extvalue == -1)
error->all(FLERR,"Fix ave/time cannot set output array intensive/extensive "
"from these inputs");
if (extarray < -1) extarray = extvalue;
else if (extvalue != extarray)
error->all(FLERR,"Fix ave/time cannot set output array intensive/extensive "
"from these inputs");
}
}
}
@ -434,32 +418,23 @@ FixAveTime::~FixAveTime()
{
// decrement lock counter in compute chunk/atom, it if still exists
if (any_variable_length &&
(nrepeat > 1 || ave == RUNNING || ave == WINDOW)) {
for (int i = 0; i < nvalues; i++)
if (varlen[i]) {
int icompute = modify->find_compute(ids[i]);
if (icompute >= 0) {
if (ave == RUNNING || ave == WINDOW)
modify->compute[icompute]->unlock(this);
modify->compute[icompute]->lock_disable();
if (any_variable_length && ((nrepeat > 1) || (ave == RUNNING) || (ave == WINDOW))) {
for (auto &val : values) {
if (val.varlen) {
auto icompute = modify->get_compute_by_id(val.id);
if (icompute) {
if ((ave == RUNNING) || (ave == WINDOW))
icompute->unlock(this);
icompute->lock_disable();
}
}
}
}
delete[] format_user;
delete[] which;
delete[] argindex;
delete[] value2index;
delete[] offcol;
delete[] varlen;
for (int i = 0; i < nvalues; i++) delete[] ids[i];
delete[] ids;
delete[] extlist;
if (fp && me == 0) {
if (fp && comm->me == 0) {
if (yaml_flag) fputs("...\n", fp);
fclose(fp);
}
@ -485,24 +460,21 @@ int FixAveTime::setmask()
void FixAveTime::init()
{
// set current indices for all computes,fixes,variables
// update indices/pointers for all computes,fixes,variables
for (int i = 0; i < nvalues; i++) {
if (which[i] == ArgInfo::COMPUTE) {
int icompute = modify->find_compute(ids[i]);
if (icompute < 0)
error->all(FLERR,"Compute ID for fix ave/time does not exist");
value2index[i] = icompute;
} else if (which[i] == ArgInfo::FIX) {
int ifix = modify->find_fix(ids[i]);
if (ifix < 0)
error->all(FLERR,"Fix ID for fix ave/time does not exist");
value2index[i] = ifix;
} else if (which[i] == ArgInfo::VARIABLE) {
int ivariable = input->variable->find(ids[i]);
if (ivariable < 0)
error->all(FLERR,"Variable name for fix ave/time does not exist");
value2index[i] = ivariable;
for (auto &val : values) {
if (val.which == ArgInfo::COMPUTE) {
val.val.c = modify->get_compute_by_id(val.id);
if (!val.val.c)
error->all(FLERR,"Compute ID {} for fix ave/time does not exist", val.id);
} else if (val.which == ArgInfo::FIX) {
val.val.f = modify->get_fix_by_id(val.id);
if (!val.val.f)
error->all(FLERR,"Fix ID {} for fix ave/time does not exist", val.id);
} else if (val.which == ArgInfo::VARIABLE) {
val.val.v = input->variable->find(val.id.c_str());
if (val.val.v < 0)
error->all(FLERR,"Variable name {} for fix ave/time does not exist", val.id);
}
}
@ -542,9 +514,6 @@ void FixAveTime::end_of_step()
void FixAveTime::invoke_scalar(bigint ntimestep)
{
int i,m;
double scalar;
// zero if first sample within single Nfreq epoch
// if any input is variable length, initialize current length
// check for exceeding length is done below
@ -556,7 +525,7 @@ void FixAveTime::invoke_scalar(bigint ntimestep)
modify->addstep_compute(ntimestep+nevery);
modify->addstep_compute(ntimestep+nfreq);
}
for (i = 0; i < nvalues; i++) vector[i] = 0.0;
for (int i = 0; i < nvalues; i++) vector[i] = 0.0;
}
// accumulate results of computes,fixes,variables to local copy
@ -564,56 +533,57 @@ void FixAveTime::invoke_scalar(bigint ntimestep)
modify->clearstep_compute();
for (i = 0; i < nvalues; i++) {
m = value2index[i];
int i = 0;
double scalar = 0.0;
for (auto &val : values) {
// invoke compute if not previously invoked
// insure no out-of-range access to variable-length compute vector
if (which[i] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[m];
if (val.which == ArgInfo::COMPUTE) {
if (argindex[i] == 0) {
if (!(compute->invoked_flag & Compute::INVOKED_SCALAR)) {
compute->compute_scalar();
compute->invoked_flag |= Compute::INVOKED_SCALAR;
if (val.argindex == 0) {
if (!(val.val.c->invoked_flag & Compute::INVOKED_SCALAR)) {
val.val.c->compute_scalar();
val.val.c->invoked_flag |= Compute::INVOKED_SCALAR;
}
scalar = compute->scalar;
scalar = val.val.c->scalar;
} else {
if (!(compute->invoked_flag & Compute::INVOKED_VECTOR)) {
compute->compute_vector();
compute->invoked_flag |= Compute::INVOKED_VECTOR;
if (!(val.val.c->invoked_flag & Compute::INVOKED_VECTOR)) {
val.val.c->compute_vector();
val.val.c->invoked_flag |= Compute::INVOKED_VECTOR;
}
if (varlen[i] && compute->size_vector < argindex[i]) scalar = 0.0;
else scalar = compute->vector[argindex[i]-1];
if (val.varlen && (val.val.c->size_vector < val.argindex)) scalar = 0.0;
else scalar = val.val.c->vector[val.argindex-1];
}
// access fix fields, guaranteed to be ready
} else if (which[i] == ArgInfo::FIX) {
if (argindex[i] == 0)
scalar = modify->fix[m]->compute_scalar();
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0)
scalar = val.val.f->compute_scalar();
else
scalar = modify->fix[m]->compute_vector(argindex[i]-1);
scalar = val.val.f->compute_vector(val.argindex-1);
// evaluate equal-style or vector-style variable
// insure no out-of-range access to vector-style variable
} else if (which[i] == ArgInfo::VARIABLE) {
if (argindex[i] == 0)
scalar = input->variable->compute_equal(m);
} else if (val.which == ArgInfo::VARIABLE) {
if (val.argindex == 0)
scalar = input->variable->compute_equal(val.val.v);
else {
double *varvec;
int nvec = input->variable->compute_vector(m,&varvec);
if (nvec < argindex[i]) scalar = 0.0;
else scalar = varvec[argindex[i]-1];
int nvec = input->variable->compute_vector(val.val.v,&varvec);
if (nvec < val.argindex) scalar = 0.0;
else scalar = varvec[val.argindex-1];
}
}
// add value to vector or just set directly if offcol is set
if (offcol[i]) vector[i] = scalar;
if (val.offcol) vector[i] = scalar;
else vector[i] += scalar;
++i;
}
// done if irepeat < nrepeat
@ -634,7 +604,7 @@ void FixAveTime::invoke_scalar(bigint ntimestep)
double repeat = nrepeat;
for (i = 0; i < nvalues; i++)
if (offcol[i] == 0) vector[i] /= repeat;
if (values[i].offcol == 0) vector[i] /= repeat;
// if ave = ONE, only single Nfreq timestep value is needed
// if ave = RUNNING, combine with all previous Nfreq timestep values
@ -667,18 +637,18 @@ void FixAveTime::invoke_scalar(bigint ntimestep)
// insure any columns with offcol set are effectively set to last value
for (i = 0; i < nvalues; i++)
if (offcol[i]) vector_total[i] = norm*vector[i];
if (values[i].offcol) vector_total[i] = norm*vector[i];
// output result to file
if (fp && me == 0) {
if (fp && comm->me == 0) {
clearerr(fp);
if (overwrite) platform::fseek(fp,filepos);
if (yaml_flag) {
if (!yaml_header || overwrite) {
yaml_header = true;
fputs("keywords: ['Step', ", fp);
for (const auto &k : keyword) fmt::print(fp, "'{}', ", k);
for (const auto &val : values) fmt::print(fp, "'{}', ", val.keyword);
fputs("]\ndata:\n", fp);
}
fmt::print(fp, " - [{}, ", ntimestep);
@ -704,8 +674,6 @@ void FixAveTime::invoke_scalar(bigint ntimestep)
void FixAveTime::invoke_vector(bigint ntimestep)
{
int i,j,m;
// first sample within single Nfreq epoch
// zero out arrays that accumulate over many samples, but not across epochs
// invoke setup_chunks() to determine current nchunk
@ -735,22 +703,20 @@ void FixAveTime::invoke_vector(bigint ntimestep)
}
int lockforever_flag = 0;
for (i = 0; i < nvalues; i++) {
if (!varlen[i] || which[i] != ArgInfo::COMPUTE) continue;
if (nrepeat > 1 && ave == ONE) {
Compute *compute = modify->compute[value2index[i]];
compute->lock(this,ntimestep,ntimestep+static_cast<bigint>(nrepeat-1)*nevery);
} else if ((ave == RUNNING || ave == WINDOW) && !lockforever) {
Compute *compute = modify->compute[value2index[i]];
compute->lock(this,update->ntimestep,-1);
for (auto &val : values) {
if (!val.varlen || (val.which != ArgInfo::COMPUTE)) continue;
if ((nrepeat > 1) && (ave == ONE)) {
val.val.c->lock(this,ntimestep,ntimestep+static_cast<bigint>(nrepeat-1)*nevery);
} else if (((ave == RUNNING) || (ave == WINDOW)) && !lockforever) {
val.val.c->lock(this,update->ntimestep,-1);
lockforever_flag = 1;
}
}
if (lockforever_flag) lockforever = 1;
}
for (i = 0; i < nrows; i++)
for (j = 0; j < nvalues; j++) array[i][j] = 0.0;
for (int i = 0; i < nrows; i++)
for (int j = 0; j < nvalues; j++) array[i][j] = 0.0;
}
// accumulate results of computes,fixes,variables to local copy
@ -758,69 +724,67 @@ void FixAveTime::invoke_vector(bigint ntimestep)
modify->clearstep_compute();
for (j = 0; j < nvalues; j++) {
m = value2index[j];
int j = 0;
for (auto &val : values) {
// invoke compute if not previously invoked
if (which[j] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[m];
if (argindex[j] == 0) {
if (!(compute->invoked_flag & Compute::INVOKED_VECTOR)) {
compute->compute_vector();
compute->invoked_flag |= Compute::INVOKED_VECTOR;
if (val.which == ArgInfo::COMPUTE) {
if (val.argindex == 0) {
if (!(val.val.c->invoked_flag & Compute::INVOKED_VECTOR)) {
val.val.c->compute_vector();
val.val.c->invoked_flag |= Compute::INVOKED_VECTOR;
}
double *cvector = compute->vector;
for (i = 0; i < nrows; i++)
double *cvector = val.val.c->vector;
for (int i = 0; i < nrows; i++)
column[i] = cvector[i];
} else {
if (!(compute->invoked_flag & Compute::INVOKED_ARRAY)) {
compute->compute_array();
compute->invoked_flag |= Compute::INVOKED_ARRAY;
if (!(val.val.c->invoked_flag & Compute::INVOKED_ARRAY)) {
val.val.c->compute_array();
val.val.c->invoked_flag |= Compute::INVOKED_ARRAY;
}
double **carray = compute->array;
int icol = argindex[j]-1;
for (i = 0; i < nrows; i++)
double **carray = val.val.c->array;
int icol = val.argindex-1;
for (int i = 0; i < nrows; i++)
column[i] = carray[i][icol];
}
// access fix fields, guaranteed to be ready
} else if (which[j] == ArgInfo::FIX) {
Fix *fix = modify->fix[m];
if (argindex[j] == 0)
for (i = 0; i < nrows; i++)
column[i] = fix->compute_vector(i);
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0)
for (int i = 0; i < nrows; i++)
column[i] = val.val.f->compute_vector(i);
else {
int icol = argindex[j]-1;
for (i = 0; i < nrows; i++)
column[i] = fix->compute_array(i,icol);
int icol = val.argindex-1;
for (int i = 0; i < nrows; i++)
column[i] = val.val.f->compute_array(i,icol);
}
// evaluate vector-style variable
// insure nvec = nrows, else error
// could be different on this timestep than when column_length(1) set nrows
} else if (which[j] == ArgInfo::VARIABLE) {
} else if (val.which == ArgInfo::VARIABLE) {
double *varvec;
int nvec = input->variable->compute_vector(m,&varvec);
int nvec = input->variable->compute_vector(val.val.v,&varvec);
if (nvec != nrows)
error->all(FLERR,"Fix ave/time vector-style variable changed length");
for (i = 0; i < nrows; i++)
error->all(FLERR,"Fix ave/time vector-style variable {} changed length", val.id);
for (int i = 0; i < nrows; i++)
column[i] = varvec[i];
}
// add columns of values to array or just set directly if offcol is set
if (offcol[j]) {
for (i = 0; i < nrows; i++)
if (val.offcol) {
for (int i = 0; i < nrows; i++)
array[i][j] = column[i];
} else {
for (i = 0; i < nrows; i++)
for (int i = 0; i < nrows; i++)
array[i][j] += column[i];
}
++j;
}
// done if irepeat < nrepeat
@ -840,38 +804,37 @@ void FixAveTime::invoke_vector(bigint ntimestep)
// unlock any variable length computes at end of Nfreq epoch
// do not unlock if ave = RUNNING or WINDOW
if (any_variable_length && nrepeat > 1 && ave == ONE) {
for (i = 0; i < nvalues; i++) {
if (!varlen[i]) continue;
Compute *compute = modify->compute[value2index[i]];
compute->unlock(this);
if (any_variable_length && (nrepeat > 1) && (ave == ONE)) {
for (auto &val : values) {
if (!val.varlen) continue;
if ((val.which == ArgInfo::COMPUTE) && val.val.c) val.val.c->unlock(this);
}
}
// average the final result for the Nfreq timestep
double repeat = nrepeat;
for (i = 0; i < nrows; i++)
for (j = 0; j < nvalues; j++)
if (offcol[j] == 0) array[i][j] /= repeat;
for (int i = 0; i < nrows; i++)
for (int j = 0; j < nvalues; j++)
if (values[j].offcol == 0) array[i][j] /= repeat;
// if ave = ONE, only single Nfreq timestep value is needed
// if ave = RUNNING, combine with all previous Nfreq timestep values
// if ave = WINDOW, combine with nwindow most recent Nfreq timestep values
if (ave == ONE) {
for (i = 0; i < nrows; i++)
for (j = 0; j < nvalues; j++) array_total[i][j] = array[i][j];
for (int i = 0; i < nrows; i++)
for (int j = 0; j < nvalues; j++) array_total[i][j] = array[i][j];
norm = 1;
} else if (ave == RUNNING) {
for (i = 0; i < nrows; i++)
for (j = 0; j < nvalues; j++) array_total[i][j] += array[i][j];
for (int i = 0; i < nrows; i++)
for (int j = 0; j < nvalues; j++) array_total[i][j] += array[i][j];
norm++;
} else if (ave == WINDOW) {
for (i = 0; i < nrows; i++)
for (j = 0; j < nvalues; j++) {
for (int i = 0; i < nrows; i++)
for (int j = 0; j < nvalues; j++) {
array_total[i][j] += array[i][j];
if (window_limit) array_total[i][j] -= array_list[iwindow][i][j];
array_list[iwindow][i][j] = array[i][j];
@ -888,32 +851,32 @@ void FixAveTime::invoke_vector(bigint ntimestep)
// insure any columns with offcol set are effectively set to last value
for (i = 0; i < nrows; i++)
for (j = 0; j < nvalues; j++)
if (offcol[j]) array_total[i][j] = norm*array[i][j];
for (int i = 0; i < nrows; i++)
for (int j = 0; j < nvalues; j++)
if (values[j].offcol) array_total[i][j] = norm*array[i][j];
// output result to file
if (fp && me == 0) {
if (fp && comm->me == 0) {
if (overwrite) platform::fseek(fp,filepos);
if (yaml_flag) {
if (!yaml_header || overwrite) {
yaml_header = true;
fputs("keywords: [", fp);
for (const auto &k : keyword) fmt::print(fp, "'{}', ", k);
for (const auto &val : values) fmt::print(fp, "'{}', ", val.keyword);
fputs("]\ndata:\n", fp);
}
fmt::print(fp, " {}:\n", ntimestep);
for (i = 0; i < nrows; i++) {
for (int i = 0; i < nrows; i++) {
fputs(" - [", fp);
for (j = 0; j < nvalues; j++) fmt::print(fp,"{}, ",array_total[i][j]/norm);
for (int j = 0; j < nvalues; j++) fmt::print(fp,"{}, ",array_total[i][j]/norm);
fputs("]\n", fp);
}
} else {
fmt::print(fp,"{} {}\n",ntimestep,nrows);
for (i = 0; i < nrows; i++) {
for (int i = 0; i < nrows; i++) {
fprintf(fp,"%d",i+1);
for (j = 0; j < nvalues; j++) fprintf(fp,format,array_total[i][j]/norm);
for (int j = 0; j < nvalues; j++) fprintf(fp,format,array_total[i][j]/norm);
fprintf(fp,"\n");
}
}
@ -938,18 +901,16 @@ int FixAveTime::column_length(int dynamic)
if (!dynamic) {
length = 0;
for (int i = 0; i < nvalues; i++) {
if (varlen[i]) continue;
if (which[i] == ArgInfo::COMPUTE) {
int icompute = modify->find_compute(ids[i]);
if (argindex[i] == 0)
lengthone = modify->compute[icompute]->size_vector;
else lengthone = modify->compute[icompute]->size_array_rows;
} else if (which[i] == ArgInfo::FIX) {
int ifix = modify->find_fix(ids[i]);
if (argindex[i] == 0) lengthone = modify->fix[ifix]->size_vector;
else lengthone = modify->fix[ifix]->size_array_rows;
} else if (which[i] == ArgInfo::VARIABLE) {
for (auto &val : values) {
if (val.varlen) continue;
if (val.which == ArgInfo::COMPUTE) {
if (val.argindex == 0)
lengthone = val.val.c->size_vector;
else lengthone = val.val.c->size_array_rows;
} else if (val.which == ArgInfo::FIX) {
if (val.argindex == 0) lengthone = val.val.f->size_vector;
else lengthone = val.val.f->size_array_rows;
} else if (val.which == ArgInfo::VARIABLE) {
// variables are always varlen = 1, so dynamic
}
if (length == 0) length = lengthone;
@ -965,15 +926,13 @@ int FixAveTime::column_length(int dynamic)
if (dynamic) {
length = 0;
for (int i = 0; i < nvalues; i++) {
if (varlen[i] == 0) continue;
m = value2index[i];
if (which[i] == ArgInfo::COMPUTE) {
Compute *compute = modify->compute[m];
lengthone = compute->lock_length();
} else if (which[i] == ArgInfo::VARIABLE) {
for (auto &val : values) {
if (val.varlen == 0) continue;
if (val.which == ArgInfo::COMPUTE) {
lengthone = val.val.c->lock_length();
} else if (val.which == ArgInfo::VARIABLE) {
double *varvec;
lengthone = input->variable->compute_vector(m,&varvec);
lengthone = input->variable->compute_vector(val.val.v,&varvec);
}
if (mode == SCALAR) continue;
if (all_variable_length) {
@ -1036,7 +995,7 @@ int FixAveTime::modify_param(int narg, char **arg)
int icol = -1;
if (utils::is_integer(arg[1])) {
icol = utils::inumeric(FLERR, arg[1], false, lmp);
if (icol < 0) icol = keyword.size() + icol + 1;
if (icol < 0) icol = values.size() + icol + 1;
icol--;
} else {
try {
@ -1045,9 +1004,9 @@ int FixAveTime::modify_param(int narg, char **arg)
icol = -1;
}
}
if ((icol < 0) || (icol >= (int) keyword.size()))
if ((icol < 0) || (icol >= (int) values.size()))
error->all(FLERR, "Thermo_modify colname column {} invalid", arg[1]);
keyword[icol] = arg[2];
values[icol].keyword = arg[2];
return 3;
}
return 0;
@ -1081,7 +1040,7 @@ void FixAveTime::options(int iarg, int narg, char **arg)
if (strcmp(arg[iarg],"file") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal fix ave/time command");
yaml_flag = utils::strmatch(arg[iarg+1],"\\.[yY][aA]?[mM][lL]$");
if (me == 0) {
if (comm->me == 0) {
fp = fopen(arg[iarg+1],"w");
if (fp == nullptr)
error->one(FLERR,"Cannot open fix ave/time file {}: {}",
@ -1140,7 +1099,7 @@ void FixAveTime::options(int iarg, int narg, char **arg)
delete[] title3;
title3 = utils::strdup(arg[iarg+1]);
iarg += 2;
} else error->all(FLERR,"Illegal fix ave/time command");
} else error->all(FLERR,"Unknown fix ave/time command option {}", arg[iarg]);
}
}

25
src/fix_ave_time.h
View File

@ -40,12 +40,24 @@ class FixAveTime : public Fix {
double compute_array(int, int) override;
private:
int me, nvalues;
int nrepeat, nfreq, irepeat;
struct value_t {
int which; // type of data: COMPUTE, FIX, VARIABLE
int argindex; // 1-based index if data is vector, else 0
int varlen; // 1 if value is from variable-length compute
int offcol;
std::string id; // compute/fix/variable ID
std::string keyword; // column keyword in output
union {
class Compute *c;
class Fix *f;
int v;
} val;
};
std::vector<value_t> values;
int nvalues, nrepeat, nfreq, irepeat;
bigint nvalid, nvalid_last;
int *which, *argindex, *value2index, *offcol;
int *varlen; // 1 if value is from variable-length compute
char **ids;
FILE *fp;
int nrows;
int any_variable_length;
@ -61,7 +73,6 @@ class FixAveTime : public Fix {
bigint filepos;
std::map<std::string, int> key2col;
std::vector<std::string> keyword;
int norm, iwindow, window_limit;
double *vector;
@ -79,8 +90,6 @@ class FixAveTime : public Fix {
void allocate_arrays();
bigint nextvalid();
};
} // namespace LAMMPS_NS
#endif
#endif