new compute global/atom command, also bug fix for descending dump sorts

This commit is contained in:
Steve Plimpton 2016-11-30 14:01:27 -07:00
parent 6759630c16
commit beb5a30f67
9 changed files with 914 additions and 5 deletions

View File

@ -767,6 +767,7 @@ KOKKOS, o = USER-OMP, t = OPT.
"erotate/sphere"_compute_erotate_sphere.html,
"erotate/sphere/atom"_compute_erotate_sphere_atom.html,
"event/displace"_compute_event_displace.html,
"global/atom"_compute_global_atom.html,
"group/group"_compute_group_group.html,
"gyration"_compute_gyration.html,
"gyration/chunk"_compute_gyration_chunk.html,

View File

@ -641,7 +641,8 @@ the restarted simulation begins.
[Related commands:]
"fix ave/chunk"_fix_ave_chunk.html
"fix ave/chunk"_fix_ave_chunk.html,
"compute global/atom"_compute_global_atom.html
[Default:]

View File

@ -0,0 +1,220 @@
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
:link(lws,http://lammps.sandia.gov)
:link(ld,Manual.html)
:link(lc,Section_commands.html#comm)
:line
compute global/atom command :h3
[Syntax:]
compute ID group-ID style index input1 input2 ... :pre
ID, group-ID are documented in "compute"_compute.html command :ulb,l
global/atom = style name of this compute command :l
index = c_ID, c_ID\[N\], f_ID, f_ID\[N\], v_name :l
c_ID = per-atom vector calculated by a compute with ID
c_ID\[I\] = Ith column of per-atom array calculated by a compute with ID
f_ID = per-atom vector calculated by a fix with ID
f_ID\[I\] = Ith column of per-atom array calculated by a fix with ID
v_name = per-atom vector calculated by an atom-style variable with name :pre
one or more inputs can be listed :l
input = c_ID, c_ID\[N\], f_ID, f_ID\[N\], v_name :l
c_ID = global vector calculated by a compute with ID
c_ID\[I\] = Ith column of global array calculated by a compute with ID, I can include wildcard (see below)
f_ID = global vector calculated by a fix with ID
f_ID\[I\] = Ith column of global array calculated by a fix with ID, I can include wildcard (see below)
v_name = global vector calculated by a vector-style variable with name :pre
:ule
[Examples:]
compute 1 all global/atom c_chunk c_com\[1\\] c_com\[2\\] c_com\[3\\]
compute 1 all global/atom c_chunk c_com\[*\\] :pre
[Description:]
Define a calculation that assigns global values to each atom from
vectors or arrays of global values. The specified {index} parameter
is used to determine which global value is assigned to each atom.
The {index} parameter must reference a per-atom vector or array from a
"compute"_compute.html or "fix"_fix.html or the evaluation of an
atom-style "variable"_variable.html. Each {input} value must
reference a global vector or array from a "compute"_compute.html or
"fix"_fix.html or the evaluation of an vector-style
"variable"_variable.html. Details are given below.
The {index} value for an atom is used as a index I (from 1 to N) into
the vector associated with each of the input values. The Ith value
from the input vector becomes one output value for that atom. If the
atom is not in the specified group, or the index I < 1 or I > M, where
M is the actual length of the input vector, then an output value of
0.0 is assigned to the atom.
An example of how this command is useful, is in the context of
"chunks" which are static or dyanmic subsets of atoms. The "compute
chunk/atom"_compute_chunk_atom.html command assigns unique chunk IDs
to each atom. It's output can be used as the {index} parameter for
this command. Various other computes with "chunk" in their style
name, such as "compute com/chunk"_compute_com_chunk.html or "compute
msd/chunk"_compute_msd_chunk.html, calculate properties for each
chunk. The output of these commands are global vectors or arrays,
with one or more values per chunk, and can be used as input values for
this command. This command will then assign the global chunk value to
each atom in the chunk, producing a per-atom vector or per-atom array
as output. The per-atom values can then be output to a dump file or
used by any command that uses per-atom values from a compute as input,
as discussed in "Section 6.15"_Section_howto.html#howto_15.
As a concrete example, these commands will calculate the displacement
of each atom from the center-of-mass of the molecule it is in, and
dump those values to a dump file. In this case, each molecule is a
chunk.
compute cc1 all chunk/atom molecule
compute myChunk all com/chunk cc1
compute prop all property/atom xu yu zu
compute glob all global/atom c_cc1 c_myChunk\[*\]
variable dx atom c_prop\[1\]-c_glob\[1\]
variable dy atom c_prop\[2\]-c_glob\[2\]
variable dz atom c_prop\[3\]-c_glob\[3\]
variable dist atom sqrt(v_dx*v_dx+v_dy*v_dy+v_dz*v_dz)
dump 1 all custom 100 tmp.dump id xu yu zu c_glob\[1\] c_glob\[2\] c_glob\[3\] &
v_dx v_dy v_dz v_dist
dump_modify 1 sort id :pre
You can add these commands to the bench/in.chain script to see how
they work.
:line
Note that for input values from a compute or fix, the bracketed index
I can be specified using a wildcard asterisk with the index to
effectively specify multiple values. This takes the form "*" or "*n"
or "n*" or "m*n". If N = the size of the vector (for {mode} = scalar)
or the number of columns in the array (for {mode} = vector), then an
asterisk with no numeric values means all indices from 1 to N. A
leading asterisk means all indices from 1 to n (inclusive). A
trailing asterisk means all indices from n to N (inclusive). A middle
asterisk means all indices from m to n (inclusive).
Using a wildcard is the same as if the individual columns of the array
had been listed one by one. E.g. these 2 compute global/atom commands
are equivalent, since the "compute com/chunk"_compute_com_chunk.html
command creates a global array with 3 columns:
compute cc1 all chunk/atom molecule
compute com all com/chunk cc1
compute 1 all global/atom c_cc1 c_com\[1\] c_com\[2\] c_com\[3\]
compute 1 all global/atom c_cc1 c_com\[*\] :pre
:line
This section explains the {index} parameter. Note that it must
reference per-atom values, as contrasted with the {input} values which
must reference global values.
Note that all of these options generate floating point values. When
they are used as an index into the specified input vectors, they
simple rounded down to convert the value to integer indices. The
final values should range from 1 to N (inclusive), since they are used
to access values from N-length vectors.
If {index} begins with "c_", a compute ID must follow which has been
previously defined in the input script. The compute must generate
per-atom quantities. See the individual "compute"_compute.html doc
page for details. If no bracketed integer is appended, the per-atom
vector calculated by the compute is used. If a bracketed integer is
appended, the Ith column of the per-atom array calculated by the
compute is used. Users can also write code for their own compute
styles and "add them to LAMMPS"_Section_modify.html. See the
discussion above for how I can be specified with a wildcard asterisk
to effectively specify multiple values.
If {index} begins with "f_", a fix ID must follow which has been
previously defined in the input script. The Fix must generate
per-atom quantities. See the individual "fix"_fix.html doc page for
details. Note that some fixes only produce their values on certain
timesteps, which must be compatible with when compute global/atom
references the values, else an error results. If no bracketed integer
is appended, the per-atom vector calculated by the fix is used. If a
bracketed integer is appended, the Ith column of the per-atom array
calculated by the fix is used. Users can also write code for their
own fix style and "add them to LAMMPS"_Section_modify.html. See the
discussion above for how I can be specified with a wildcard asterisk
to effectively specify multiple values.
If {index} begins with "v_", a variable name must follow which has
been previously defined in the input script. It must be an
"atom-style variable"_variable.html. Atom-style variables can
reference thermodynamic keywords and various per-atom attributes, or
invoke other computes, fixes, or variables when they are evaluated, so
this is a very general means of generating per-atom quantities to use
as {index}.
:line
This section explains the kinds of {input} values that can be used.
Note that inputs reference global values, as contrasted with the
{index} parameter which must reference per-atom values.
If a value begins with "c_", a compute ID must follow which has been
previously defined in the input script. The compute must generate a
global vector or array. See the individual "compute"_compute.html doc
page for details. If no bracketed integer is appended, the vector
calculated by the compute is used. If a bracketed integer is
appended, the Ith column of the array calculated by the compute is
used. Users can also write code for their own compute styles and "add
them to LAMMPS"_Section_modify.html. See the discussion above for how
I can be specified with a wildcard asterisk to effectively specify
multiple values.
If a value begins with "f_", a fix ID must follow which has been
previously defined in the input script. The fix must generate a
global vector or array. See the individual "fix"_fix.html doc page
for details. Note that some fixes only produce their values on
certain timesteps, which must be compatible with when compute
global/atom references the values, else an error results. If no
bracketed integer is appended, the vector calculated by the fix is
used. If a bracketed integer is appended, the Ith column of the array
calculated by the fix is used. Users can also write code for their
own fix style and "add them to LAMMPS"_Section_modify.html. See the
discussion above for how I can be specified with a wildcard asterisk
to effectively specify multiple values.
If a value begins with "v_", a variable name must follow which has
been previously defined in the input script. It must be a
"vector-style variable"_variable.html. Vector-style variables can
reference thermodynamic keywords and various other attributes of
atoms, or invoke other computes, fixes, or variables when they are
evaluated, so this is a very general means of generating a vector of
global quantities which the {index} parameter will reference for
assignement of global values to atoms.
:line
[Output info:]
If a single input is specified this compute produces a per-atom
vector. If multiple inputs are specified, this compute produces a
per-atom array values, where the number of columns is equal to the
number of inputs specified. These values can be used by any command
that uses per-atom vector or array values from a compute as input.
See "Section 6.15"_Section_howto.html#howto_15 for an overview of
LAMMPS output options.
The per-atom vector or array values will be in whatever units the
corresponsing input values are in.
[Restrictions:] none
[Related commands:]
"compute"_compute.html, "fix"_fix.html, "variable"_variable.html,
"compute chunk/atom"_compute_chunk_atom.html, "compute
reduce"_compute_reduce.html
[Default:] none

View File

@ -311,6 +311,7 @@ compute_erotate_sphere.html
compute_erotate_sphere_atom.html
compute_event_displace.html
compute_fep.html
compute_global_atom.html
compute_group_group.html
compute_gyration.html
compute_gyration_chunk.html

View File

@ -40,7 +40,8 @@ int Compute::instance_total = 0;
Compute::Compute(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp),
id(NULL), style(NULL),
vector(NULL), array(NULL), vector_atom(NULL), array_atom(NULL), vector_local(NULL), array_local(NULL),
vector(NULL), array(NULL), vector_atom(NULL),
array_atom(NULL), vector_local(NULL), array_local(NULL),
tlist(NULL), vbiasall(NULL)
{
instance_me = instance_total++;

View File

@ -70,7 +70,7 @@ ComputeBondLocal::ComputeBondLocal(LAMMPS *lmp, int narg, char **arg) :
// set velflag if compute any quantities based on velocities
singleflag = 0;
ghostvelflag = 0;
velflag = 0;
for (int i = 0; i < nvalues; i++) {
if (bstyle[i] == ENGPOT || bstyle[i] == FORCE) singleflag = 1;
if (bstyle[i] == VELVIB || bstyle[i] == OMEGA || bstyle[i] == ENGTRANS ||

540
src/compute_global_atom.cpp Normal file
View File

@ -0,0 +1,540 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#include <string.h>
#include <stdlib.h>
#include "compute_global_atom.h"
#include "atom.h"
#include "update.h"
#include "domain.h"
#include "modify.h"
#include "fix.h"
#include "force.h"
#include "comm.h"
#include "group.h"
#include "input.h"
#include "variable.h"
#include "memory.h"
#include "error.h"
using namespace LAMMPS_NS;
enum{COMPUTE,FIX,VARIABLE};
enum{VECTOR,ARRAY};
#define INVOKED_VECTOR 2
#define INVOKED_ARRAY 4
#define INVOKED_PERATOM 8
#define BIG 1.0e20
/* ---------------------------------------------------------------------- */
ComputeGlobalAtom::ComputeGlobalAtom(LAMMPS *lmp, int narg, char **arg) :
Compute(lmp, narg, arg),
idref(NULL), which(NULL), argindex(NULL), value2index(NULL), ids(NULL),
indices(NULL), varatom(NULL), vecglobal(NULL)
{
if (narg < 5) error->all(FLERR,"Illegal compute global/atom command");
peratom_flag = 1;
// process index arg
int iarg = 3;
if (strncmp(arg[iarg],"c_",2) == 0 ||
strncmp(arg[iarg],"f_",2) == 0 ||
strncmp(arg[iarg],"v_",2) == 0) {
if (arg[iarg][0] == 'c') whichref = COMPUTE;
else if (arg[iarg][0] == 'f') whichref = FIX;
else if (arg[iarg][0] == 'v') whichref = VARIABLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
strcpy(suffix,&arg[iarg][2]);
char *ptr = strchr(suffix,'[');
if (ptr) {
if (suffix[strlen(suffix)-1] != ']')
error->all(FLERR,"Illegal compute global/atom command");
indexref = atoi(ptr+1);
*ptr = '\0';
} else indexref = 0;
n = strlen(suffix) + 1;
idref = new char[n];
strcpy(idref,suffix);
delete [] suffix;
} else error->all(FLERR,"Illegal compute global/atom command");
iarg++;
// expand args if any have wildcard character "*"
int expand = 0;
char **earg;
int nargnew = input->expand_args(narg-iarg,&arg[iarg],1,earg);
if (earg != &arg[iarg]) expand = 1;
arg = earg;
// parse values until one isn't recognized
which = new int[nargnew];
argindex = new int[nargnew];
ids = new char*[nargnew];
value2index = new int[nargnew];
nvalues = 0;
iarg = 0;
while (iarg < nargnew) {
ids[nvalues] = NULL;
if (strncmp(arg[iarg],"c_",2) == 0 ||
strncmp(arg[iarg],"f_",2) == 0 ||
strncmp(arg[iarg],"v_",2) == 0) {
if (arg[iarg][0] == 'c') which[nvalues] = COMPUTE;
else if (arg[iarg][0] == 'f') which[nvalues] = FIX;
else if (arg[iarg][0] == 'v') which[nvalues] = VARIABLE;
int n = strlen(arg[iarg]);
char *suffix = new char[n];
strcpy(suffix,&arg[iarg][2]);
char *ptr = strchr(suffix,'[');
if (ptr) {
if (suffix[strlen(suffix)-1] != ']')
error->all(FLERR,"Illegal compute global/atom command");
argindex[nvalues] = atoi(ptr+1);
*ptr = '\0';
} else argindex[nvalues] = 0;
n = strlen(suffix) + 1;
ids[nvalues] = new char[n];
strcpy(ids[nvalues],suffix);
nvalues++;
delete [] suffix;
} else error->all(FLERR,"Illegal compute global/atom command");
iarg++;
}
// if wildcard expansion occurred, free earg memory from expand_args()
if (expand) {
for (int i = 0; i < nargnew; i++) delete [] earg[i];
memory->sfree(earg);
}
// setup and error check both index arg and values
if (whichref == COMPUTE) {
int icompute = modify->find_compute(idref);
if (icompute < 0)
error->all(FLERR,"Compute ID for compute global/atom does not exist");
if (!modify->compute[icompute]->peratom_flag)
error->all(FLERR,"Compute global/atom compute does not "
"calculate a per-atom vector or array");
if (indexref == 0 &&
modify->compute[icompute]->size_peratom_cols != 0)
error->all(FLERR,"Compute global/atom compute does not "
"calculate a per-atom vector");
if (indexref && modify->compute[icompute]->size_peratom_cols == 0)
error->all(FLERR,"Compute global/atom compute does not "
"calculate a per-atom array");
if (indexref && indexref > modify->compute[icompute]->size_peratom_cols)
error->all(FLERR,
"Compute global/atom compute array is accessed out-of-range");
} else if (whichref == FIX) {
int ifix = modify->find_fix(idref);
if (ifix < 0)
error->all(FLERR,"Fix ID for compute global/atom does not exist");
if (!modify->fix[ifix]->peratom_flag)
error->all(FLERR,"Compute global/atom fix does not "
"calculate a per-atom vector or array");
if (indexref == 0 &&
modify->fix[ifix]->size_peratom_cols != 0)
error->all(FLERR,"Compute global/atom fix does not "
"calculate a per-atom vector");
if (indexref && modify->fix[ifix]->size_peratom_cols == 0)
error->all(FLERR,"Compute global/atom fix does not "
"calculate a per-atom array");
if (indexref && indexref > modify->fix[ifix]->size_peratom_cols)
error->all(FLERR,
"Compute global/atom fix array is accessed out-of-range");
} else if (whichref == VARIABLE) {
int ivariable = input->variable->find(idref);
if (ivariable < 0)
error->all(FLERR,"Variable name for compute global/atom does not exist");
if (input->variable->atomstyle(ivariable) == 0)
error->all(FLERR,"Compute global/atom variable is not "
"atom-style variable");
}
for (int i = 0; i < nvalues; i++) {
if (which[i] == COMPUTE) {
int icompute = modify->find_compute(ids[i]);
if (icompute < 0)
error->all(FLERR,"Compute ID for compute global/atom does not exist");
if (argindex[i] == 0) {
if (!modify->compute[icompute]->vector_flag)
error->all(FLERR,"Compute global/atom compute does not "
"calculate a global vector");
} else {
if (!modify->compute[icompute]->array_flag)
error->all(FLERR,"Compute global/atom compute does not "
"calculate a global array");
if (argindex[i] > modify->compute[icompute]->size_array_cols)
error->all(FLERR,"Compute global/atom compute array is "
"accessed out-of-range");
}
} else if (which[i] == FIX) {
int ifix = modify->find_fix(ids[i]);
if (ifix < 0)
error->all(FLERR,"Fix ID for compute global/atom does not exist");
if (argindex[i] == 0) {
if (!modify->fix[ifix]->vector_flag)
error->all(FLERR,"Compute global/atom fix does not "
"calculate a global vector");
} else {
if (!modify->fix[ifix]->array_flag)
error->all(FLERR,"Compute global/atom fix does not "
"calculate a global array");
if (argindex[i] > modify->fix[ifix]->size_array_cols)
error->all(FLERR,"Compute global/atom fix array is "
"accessed out-of-range");
}
} else if (which[i] == VARIABLE) {
int ivariable = input->variable->find(ids[i]);
if (ivariable < 0)
error->all(FLERR,"Variable name for compute global/atom "
"does not exist");
if (input->variable->vectorstyle(ivariable) == 0)
error->all(FLERR,"Compute global/atom variable is not "
"vector-style variable");
}
}
// this compute produces either a peratom vector or array
if (nvalues == 1) size_peratom_cols = 0;
else size_peratom_cols = nvalues;
nmax = maxvector = 0;
indices = NULL;
varatom = NULL;
vecglobal = NULL;
vector_atom = NULL;
array_atom = NULL;
}
/* ---------------------------------------------------------------------- */
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(varatom);
memory->destroy(vecglobal);
memory->destroy(vector_atom);
memory->destroy(array_atom);
}
/* ---------------------------------------------------------------------- */
void ComputeGlobalAtom::init()
{
// set indices of all computes,fixes,variables
if (whichref == COMPUTE) {
int icompute = modify->find_compute(idref);
if (icompute < 0)
error->all(FLERR,"Compute ID for compute global/atom does not exist");
ref2index = icompute;
} else if (whichref == FIX) {
int ifix = modify->find_fix(idref);
if (ifix < 0)
error->all(FLERR,"Fix ID for compute global/atom does not exist");
ref2index = ifix;
} else if (whichref == VARIABLE) {
int ivariable = input->variable->find(idref);
if (ivariable < 0)
error->all(FLERR,"Variable name for compute global/atom does not exist");
ref2index = ivariable;
}
for (int m = 0; m < nvalues; m++) {
if (which[m] == COMPUTE) {
int icompute = modify->find_compute(ids[m]);
if (icompute < 0)
error->all(FLERR,"Compute ID for compute global/atom does not exist");
value2index[m] = icompute;
} else if (which[m] == FIX) {
int ifix = modify->find_fix(ids[m]);
if (ifix < 0)
error->all(FLERR,"Fix ID for compute global/atom does not exist");
value2index[m] = ifix;
} else if (which[m] == VARIABLE) {
int ivariable = input->variable->find(ids[m]);
if (ivariable < 0)
error->all(FLERR,"Variable name for compute global/atom "
"does not exist");
value2index[m] = ivariable;
}
}
}
/* ---------------------------------------------------------------------- */
void ComputeGlobalAtom::compute_peratom()
{
int i,j;
invoked_peratom = update->ntimestep;
// grow indices and output vector or array if necessary
if (atom->nmax > nmax) {
nmax = atom->nmax;
memory->destroy(indices);
memory->create(indices,nmax,"global/atom:indices");
if (nvalues == 1) {
memory->destroy(vector_atom);
memory->create(vector_atom,nmax,"global/atom:vector_atom");
} else {
memory->destroy(array_atom);
memory->create(array_atom,nmax,nvalues,"global/atom:array_atom");
}
}
// setup current peratom indices
// integer indices are rounded down from double values
// indices are decremented from 1 to N -> 0 to N-1
int *mask = atom->mask;
int nlocal = atom->nlocal;
if (whichref == COMPUTE) {
Compute *compute = modify->compute[ref2index];
if (!(compute->invoked_flag & INVOKED_PERATOM)) {
compute->compute_peratom();
compute->invoked_flag |= INVOKED_PERATOM;
}
if (indexref == 0) {
double *compute_vector = compute->vector_atom;
for (i = 0; i < nlocal; i++)
if (mask[i] & groupbit)
indices[i] = static_cast<int> (compute_vector[i]) - 1;
} else {
double **compute_array = compute->array_atom;
int im1 = indexref - 1;
for (i = 0; i < nlocal; i++)
if (mask[i] & groupbit)
indices[i] = static_cast<int> (compute_array[i][im1]) - 1;
}
} else if (whichref == FIX) {
if (update->ntimestep % modify->fix[ref2index]->peratom_freq)
error->all(FLERR,"Fix used in compute global/atom not "
"computed at compatible time");
Fix *fix = modify->fix[ref2index];
if (indexref == 0) {
double *fix_vector = fix->vector_atom;
for (i = 0; i < nlocal; i++)
if (mask[i] & groupbit)
indices[i] = static_cast<int> (fix_vector[i]) - 1;
} else {
double **fix_array = fix->array_atom;
int im1 = indexref - 1;
for (i = 0; i < nlocal; i++)
if (mask[i] & groupbit)
indices[i] = static_cast<int> (fix_array[i][im1]) - 1;
}
} else if (whichref == VARIABLE) {
if (atom->nmax > nmax) {
nmax = atom->nmax;
memory->destroy(varatom);
memory->create(varatom,nmax,"global/atom:varatom");
}
input->variable->compute_atom(ref2index,igroup,varatom,1,0);
for (i = 0; i < nlocal; i++)
if (mask[i] & groupbit)
indices[i] = static_cast<int> (varatom[i]) - 1;
}
// loop over values to fill output vector or array
for (int m = 0; m < nvalues; m++) {
// output = vector
if (argindex[m] == 0) {
int vmax;
double *source;
if (which[m] == COMPUTE) {
Compute *compute = modify->compute[value2index[m]];
if (!(compute->invoked_flag & INVOKED_VECTOR)) {
compute->compute_vector();
compute->invoked_flag |= INVOKED_VECTOR;
}
source = compute->vector;
vmax = compute->size_vector;
} else if (which[m] == FIX) {
if (update->ntimestep % modify->fix[value2index[m]]->peratom_freq)
error->all(FLERR,"Fix used in compute global/atom not "
"computed at compatible time");
Fix *fix = modify->fix[value2index[m]];
vmax = fix->size_vector;
if (vmax > maxvector) {
maxvector = vmax;
memory->destroy(vecglobal);
memory->create(vecglobal,maxvector,"global/atom:vecglobal");
}
for (i = 0; i < vmax; i++)
vecglobal[i] = fix->compute_vector(i);
source = vecglobal;
} else if (which[m] == VARIABLE) {
vmax = input->variable->compute_vector(value2index[m],&source);
}
if (nvalues == 1) {
for (i = 0; i < nlocal; i++) {
vector_atom[i] = 0.0;
if (mask[i] & groupbit) {
j = indices[i];
if (j >= 0 && j < vmax) vector_atom[i] = source[j];
}
}
} else {
for (i = 0; i < nlocal; i++) {
array_atom[i][m] = 0.0;
if (mask[i] & groupbit) {
j = indices[i];
if (j >= 0 && j < vmax) array_atom[i][m] = source[j];
}
}
}
// output = array
} else {
int vmax;
double *source;
int col = argindex[m] - 1;
if (which[m] == COMPUTE) {
Compute *compute = modify->compute[value2index[m]];
if (!(compute->invoked_flag & INVOKED_ARRAY)) {
compute->compute_array();
compute->invoked_flag |= INVOKED_ARRAY;
}
double **compute_array = compute->array;
vmax = compute->size_array_rows;
if (vmax > maxvector) {
maxvector = vmax;
memory->destroy(vecglobal);
memory->create(vecglobal,maxvector,"global/atom:vecglobal");
}
for (i = 0; i < vmax; i++)
vecglobal[i] = compute_array[i][col];
source = vecglobal;
} else if (which[m] == FIX) {
if (update->ntimestep % modify->fix[value2index[m]]->peratom_freq)
error->all(FLERR,"Fix used in compute global/atom not "
"computed at compatible time");
Fix *fix = modify->fix[value2index[m]];
vmax = fix->size_array_rows;
if (vmax > maxvector) {
maxvector = vmax;
memory->destroy(vecglobal);
memory->create(vecglobal,maxvector,"global/atom:vecglobal");
}
for (i = 0; i < vmax; i++)
vecglobal[i] = fix->compute_array(i,col);
source = vecglobal;
} else if (which[m] == VARIABLE) {
vmax = input->variable->compute_vector(value2index[m],&source);
}
if (nvalues == 1) {
for (i = 0; i < nlocal; i++) {
vector_atom[i] = 0.0;
if (mask[i] & groupbit) {
j = indices[i];
if (j >= 0 && j < vmax) vector_atom[i] = source[j];
}
}
} else {
for (i = 0; i < nlocal; i++) {
array_atom[i][m] = 0.0;
if (mask[i] & groupbit) {
j = indices[i];
if (j >= 0 && j < vmax) array_atom[i][m] = source[j];
}
}
}
}
}
}
/* ----------------------------------------------------------------------
memory usage of local atom-based array
------------------------------------------------------------------------- */
double ComputeGlobalAtom::memory_usage()
{
double bytes = nmax*nvalues * sizeof(double);
bytes += nmax * sizeof(int); // indices
if (varatom) bytes += nmax * sizeof(double); // varatom
bytes += maxvector * sizeof(double); // vecglobal
return bytes;
}

143
src/compute_global_atom.h Normal file
View File

@ -0,0 +1,143 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifdef COMPUTE_CLASS
ComputeStyle(global/atom,ComputeGlobalAtom)
#else
#ifndef LMP_COMPUTE_GLOBAL_ATOM_H
#define LMP_COMPUTE_GLOBAL_ATOM_H
#include "compute.h"
namespace LAMMPS_NS {
class ComputeGlobalAtom : public Compute {
public:
ComputeGlobalAtom(class LAMMPS *, int, char **);
virtual ~ComputeGlobalAtom();
void init();
void compute_peratom();
double memory_usage();
protected:
int whichref,indexref,ref2index;
char *idref;
int nvalues;
int *which,*argindex,*value2index;
char **ids;
int nmax,maxvector;
int *indices;
double *varatom;
double *vecglobal;
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Illegal ... command
Self-explanatory. Check the input script syntax and compare to the
documentation for the command. You can use -echo screen as a
command-line option when running LAMMPS to see the offending line.
E: Region ID for compute reduce/region does not exist
Self-explanatory.
E: Compute reduce replace requires min or max mode
Self-explanatory.
E: Invalid replace values in compute reduce
Self-explanatory.
E: Compute ID for compute reduce does not exist
Self-explanatory.
E: Compute reduce compute does not calculate a per-atom vector
Self-explanatory.
E: Compute reduce compute does not calculate a per-atom array
Self-explanatory.
E: Compute reduce compute array is accessed out-of-range
An index for the array is out of bounds.
E: Compute reduce compute does not calculate a local vector
Self-explanatory.
E: Compute reduce compute does not calculate a local array
Self-explanatory.
E: Compute reduce compute calculates global values
A compute that calculates peratom or local values is required.
E: Fix ID for compute reduce does not exist
Self-explanatory.
E: Compute reduce fix does not calculate a per-atom vector
Self-explanatory.
E: Compute reduce fix does not calculate a per-atom array
Self-explanatory.
E: Compute reduce fix array is accessed out-of-range
An index for the array is out of bounds.
E: Compute reduce fix does not calculate a local vector
Self-explanatory.
E: Compute reduce fix does not calculate a local array
Self-explanatory.
E: Compute reduce fix calculates global values
A fix that calculates peratom or local values is required.
E: Variable name for compute reduce does not exist
Self-explanatory.
E: Compute reduce variable is not atom-style variable
Self-explanatory.
E: Fix used in compute reduce not computed at compatible time
Fixes generate their values on specific timesteps. Compute reduce is
requesting a value on a non-allowed timestep.
*/

View File

@ -43,8 +43,10 @@ enum{SCALAR,VECTOR};
FixAveTime::FixAveTime(LAMMPS *lmp, int narg, char **arg) :
Fix(lmp, narg, arg),
nvalues(0), which(NULL), argindex(NULL), value2index(NULL), offcol(NULL), varlen(NULL), ids(NULL),
fp(NULL), offlist(NULL), format(NULL), format_user(NULL), vector(NULL), vector_total(NULL), vector_list(NULL),
nvalues(0), which(NULL), argindex(NULL), value2index(NULL),
offcol(NULL), varlen(NULL), ids(NULL),
fp(NULL), offlist(NULL), format(NULL), format_user(NULL),
vector(NULL), vector_total(NULL), vector_list(NULL),
column(NULL), array(NULL), array_total(NULL), array_list(NULL)
{
if (narg < 7) error->all(FLERR,"Illegal fix ave/time command");