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

Merge pull request #3052 from lammps/time-dumps2

Browse Source 
Time-dependent dumps for variable timestep (alternate implementation)
This commit is contained in:
Axel Kohlmeyer
2021-12-23 16:34:12 -05:00
committed by GitHub
parent 6af36075ba a653ee6b2c
commit b39d1993bb
18 changed files with 857 additions and 436 deletions
Download Patch File Download Diff File Expand all files Collapse all files

78
doc/src/dump.rst
View File

@ -169,11 +169,12 @@ or multiple smaller files).
.. note:: .. note::
Because periodic boundary conditions are enforced only on Because periodic boundary conditions are enforced only on timesteps
timesteps when neighbor lists are rebuilt, the coordinates of an atom when neighbor lists are rebuilt, the coordinates of an atom written
written to a dump file may be slightly outside the simulation box. to a dump file may be slightly outside the simulation box.
Re-neighbor timesteps will not typically coincide with the timesteps Re-neighbor timesteps will not typically coincide with the
dump snapshots are written. See the :doc:`dump_modify pbc <dump_modify>` command if you with to force coordinates to be timesteps dump snapshots are written. See the :doc:`dump_modify
pbc <dump_modify>` command if you with to force coordinates to be
strictly inside the simulation box. strictly inside the simulation box.
.. note:: .. note::
@ -189,20 +190,21 @@ or multiple smaller files).
multiple processors, each of which owns a subset of the atoms. multiple processors, each of which owns a subset of the atoms.
For the *atom*, *custom*, *cfg*, and *local* styles, sorting is off by For the *atom*, *custom*, *cfg*, and *local* styles, sorting is off by
default. For the *dcd*, *xtc*, *xyz*, and *molfile* styles, sorting by default. For the *dcd*, *xtc*, *xyz*, and *molfile* styles, sorting
atom ID is on by default. See the :doc:`dump_modify <dump_modify>` doc by atom ID is on by default. See the :doc:`dump_modify <dump_modify>`
page for details. doc page for details.
The *atom/gz*, *cfg/gz*, *custom/gz*, *local/gz*, and *xyz/gz* styles are identical The *atom/gz*, *cfg/gz*, *custom/gz*, *local/gz*, and *xyz/gz* styles
in command syntax to the corresponding styles without "gz", however, are identical in command syntax to the corresponding styles without
they generate compressed files using the zlib library. Thus the filename "gz", however, they generate compressed files using the zlib
suffix ".gz" is mandatory. This is an alternative approach to writing library. Thus the filename suffix ".gz" is mandatory. This is an
compressed files via a pipe, as done by the regular dump styles, which alternative approach to writing compressed files via a pipe, as done
may be required on clusters where the interface to the high-speed network by the regular dump styles, which may be required on clusters where
disallows using the fork() library call (which is needed for a pipe). the interface to the high-speed network disallows using the fork()
For the remainder of this doc page, you should thus consider the *atom* library call (which is needed for a pipe). For the remainder of this
and *atom/gz* styles (etc) to be inter-changeable, with the exception doc page, you should thus consider the *atom* and *atom/gz* styles
of the required filename suffix. (etc) to be inter-changeable, with the exception of the required
filename suffix.
Similarly, the *atom/zstd*, *cfg/zstd*, *custom/zstd*, *local/zstd*, Similarly, the *atom/zstd*, *cfg/zstd*, *custom/zstd*, *local/zstd*,
and *xyz/zstd* styles are identical to the gz styles, but use the Zstd and *xyz/zstd* styles are identical to the gz styles, but use the Zstd
@ -280,10 +282,11 @@ This bounding box is convenient for many visualization programs. The
meaning of the 6 character flags for "xx yy zz" is the same as above. meaning of the 6 character flags for "xx yy zz" is the same as above.
Note that the first two numbers on each line are now xlo_bound instead Note that the first two numbers on each line are now xlo_bound instead
of xlo, etc, since they represent a bounding box. See the :doc:`Howto triclinic <Howto_triclinic>` page for a geometric description of xlo, etc, since they represent a bounding box. See the :doc:`Howto
of triclinic boxes, as defined by LAMMPS, simple formulas for how the triclinic <Howto_triclinic>` page for a geometric description of
6 bounding box extents (xlo_bound,xhi_bound,etc) are calculated from triclinic boxes, as defined by LAMMPS, simple formulas for how the 6
the triclinic parameters, and how to transform those parameters to and bounding box extents (xlo_bound,xhi_bound,etc) are calculated from the
triclinic parameters, and how to transform those parameters to and
from other commonly used triclinic representations. from other commonly used triclinic representations.
The "ITEM: ATOMS" line in each snapshot lists column descriptors for The "ITEM: ATOMS" line in each snapshot lists column descriptors for
@ -315,23 +318,24 @@ written to the dump file. This local data is typically calculated by
each processor based on the atoms it owns, but there may be zero or each processor based on the atoms it owns, but there may be zero or
more entities per atom, e.g. a list of bond distances. An explanation more entities per atom, e.g. a list of bond distances. An explanation
of the possible dump local attributes is given below. Note that by of the possible dump local attributes is given below. Note that by
using input from the :doc:`compute property/local <compute_property_local>` command with dump local, using input from the :doc:`compute property/local
it is possible to generate information on bonds, angles, etc that can <compute_property_local>` command with dump local, it is possible to
be cut and pasted directly into a data file read by the generate information on bonds, angles, etc that can be cut and pasted
:doc:`read_data <read_data>` command. directly into a data file read by the :doc:`read_data <read_data>`
command.
Style *cfg* has the same command syntax as style *custom* and writes Style *cfg* has the same command syntax as style *custom* and writes
extended CFG format files, as used by the extended CFG format files, as used by the `AtomEye
`AtomEye <http://li.mit.edu/Archive/Graphics/A/>`_ visualization <http://li.mit.edu/Archive/Graphics/A/>`_ visualization package.
package. Since the extended CFG format uses a single snapshot of the Since the extended CFG format uses a single snapshot of the system per
system per file, a wildcard "\*" must be included in the filename, as file, a wildcard "\*" must be included in the filename, as discussed
discussed below. The list of atom attributes for style *cfg* must below. The list of atom attributes for style *cfg* must begin with
begin with either "mass type xs ys zs" or "mass type xsu ysu zsu" either "mass type xs ys zs" or "mass type xsu ysu zsu" since these
since these quantities are needed to write the CFG files in the quantities are needed to write the CFG files in the appropriate format
appropriate format (though the "mass" and "type" fields do not appear (though the "mass" and "type" fields do not appear explicitly in the
explicitly in the file). Any remaining attributes will be stored as file). Any remaining attributes will be stored as "auxiliary
"auxiliary properties" in the CFG files. Note that you will typically properties" in the CFG files. Note that you will typically want to
want to use the :doc:`dump_modify element <dump_modify>` command with use the :doc:`dump_modify element <dump_modify>` command with
CFG-formatted files, to associate element names with atom types, so CFG-formatted files, to associate element names with atom types, so
that AtomEye can render atoms appropriately. When unwrapped that AtomEye can render atoms appropriately. When unwrapped
coordinates *xsu*, *ysu*, and *zsu* are requested, the nominal AtomEye coordinates *xsu*, *ysu*, and *zsu* are requested, the nominal AtomEye

151
doc/src/dump_modify.rst
View File

@ -17,7 +17,7 @@ Syntax
* one or more keyword/value pairs may be appended * one or more keyword/value pairs may be appended
* these keywords apply to various dump styles * these keywords apply to various dump styles
* keyword = *append* or *at* or *buffer* or *delay* or *element* or *every* or *fileper* or *first* or *flush* or *format* or *header* or *image* or *label* or *maxfiles* or *nfile* or *pad* or *pbc* or *precision* or *region* or *refresh* or *scale* or *sfactor* or *sort* or *tfactor* or *thermo* or *thresh* or *time* or *units* or *unwrap* * keyword = *append* or *at* or *buffer* or *delay* or *element* or *every* or *every/time* or *fileper* or *first* or *flush* or *format* or *header* or *image* or *label* or *maxfiles* or *nfile* or *pad* or *pbc* or *precision* or *region* or *refresh* or *scale* or *sfactor* or *sort* or *tfactor* or *thermo* or *thresh* or *time* or *units* or *unwrap*
.. parsed-literal:: .. parsed-literal::
@ -32,6 +32,9 @@ Syntax
*every* arg = N *every* arg = N
N = dump every this many timesteps N = dump every this many timesteps
N can be a variable (see below) N can be a variable (see below)
*every/time* arg = Delta
Delta = dump every this interval in simulation time (time units)
Delta can be a variable (see below)
*fileper* arg = Np *fileper* arg = Np
Np = write one file for every this many processors Np = write one file for every this many processors
*first* arg = *yes* or *no* *first* arg = *yes* or *no*
@ -197,11 +200,19 @@ will be accepted.
---------- ----------
The *every* keyword changes the dump frequency originally specified by The *every* keyword can be used with any dump style except the *dcd*
the :doc:`dump <dump>` command to a new value. The every keyword can be and *xtc* styles. It does two things. It specifies that the interval
specified in one of two ways. It can be a numeric value in which case between dump snapshots will be set in timesteps, which is the default
it must be > 0. Or it can be an :doc:`equal-style variable <variable>`, if the *every* or *every/time* keywords are not used. See the
which should be specified as v_name, where name is the variable name. *every/time* keyword for how to specify the interval in simulation
time, i.e. in time units of the :doc:`units <units>` command. The
*every* keyword also sets the interval value, which overrides the dump
frequency originally specified by the :doc:`dump <dump>` command.
The *every* keyword can be specified in one of two ways. It can be a
numeric value in which case it must be > 0. Or it can be an
:doc:`equal-style variable <variable>`, which should be specified as
v_name, where name is the variable name.
In this case, the variable is evaluated at the beginning of a run to In this case, the variable is evaluated at the beginning of a run to
determine the next timestep at which a dump snapshot will be written determine the next timestep at which a dump snapshot will be written
@ -210,11 +221,12 @@ determine the next timestep, etc. Thus the variable should return
timestep values. See the stagger() and logfreq() and stride() math timestep values. See the stagger() and logfreq() and stride() math
functions for :doc:`equal-style variables <variable>`, as examples of functions for :doc:`equal-style variables <variable>`, as examples of
useful functions to use in this context. Other similar math functions useful functions to use in this context. Other similar math functions
could easily be added as options for :doc:`equal-style variables <variable>`. Also see the next() function, which allows could easily be added as options for :doc:`equal-style variables
use of a file-style variable which reads successive values from a <variable>`. Also see the next() function, which allows use of a
file, each time the variable is evaluated. Used with the *every* file-style variable which reads successive values from a file, each
keyword, if the file contains a list of ascending timesteps, you can time the variable is evaluated. Used with the *every* keyword, if the
output snapshots whenever you wish. file contains a list of ascending timesteps, you can output snapshots
whenever you wish.
Note that when using the variable option with the *every* keyword, you Note that when using the variable option with the *every* keyword, you
need to use the *first* option if you want an initial snapshot written need to use the *first* option if you want an initial snapshot written
@ -255,14 +267,103 @@ in file tmp.times:
---------- ----------
The *every/time* keyword can be used with any dump style except the
*dcd* and *xtc* styles. It does two things. It specifies that the
interval between dump snapshots will be set in simulation time,
i.e. in time units of the :doc:`units <units>` command. This can be
useful when the timestep size varies during a simulation run, e.g. by
use of the :doc:`fix dt/reset <fix_dt_reset>` command. The default is
to specify the interval in timesteps; see the *every* keyword. The
*every/time* command also sets the interval value.
.. note::
If you wish dump styles *atom*, *custom*, *local*, or *xyz* to
include the simulation time as a field in the header portion of
each snapshot, you also need to use the dump_modify *time* keyword
with a setting of *yes*. See its documentation below.
Note that since snapshots are output on simulation steps, each
snapshot will be written on the first timestep whose associated
simulation time is >= the exact snapshot time value.
As with the *every* option, the *Delta* value can be specified in one
of two ways. It can be a numeric value in which case it must be >
0.0. Or it can be an :doc:`equal-style variable <variable>`, which
should be specified as v_name, where name is the variable name.
In this case, the variable is evaluated at the beginning of a run to
determine the next simulation time at which a dump snapshot will be
written out. On that timestep the variable will be evaluated again to
determine the next simulation time, etc. Thus the variable should
return values in time units. Note the current timestep or simulation
time can be used in an :doc:`equal-style variables <variable>` since
they are both thermodynamic keywords. Also see the next() function,
which allows use of a file-style variable which reads successive
values from a file, each time the variable is evaluated. Used with
the *every/time* keyword, if the file contains a list of ascending
simulation times, you can output snapshots whenever you wish.
Note that when using the variable option with the *every/time*
keyword, you need to use the *first* option if you want an initial
snapshot written to the dump file. The *every/time* keyword cannot be
used with the dump *dcd* style.
For example, the following commands will write snapshots at successive
simulation times which grow by a factor of 1.5 with each interval.
The dt value used in the variable is to avoid a zero result when the
initial simulation time is 0.0.
.. code-block:: LAMMPS
variable increase equal 1.5*(time+dt)
dump 1 all atom 100 tmp.dump
dump_modify 1 every/time v_increase first yes
The following commands would write snapshots at the times listed in
file tmp.times:
.. code-block:: LAMMPS
variable f file tmp.times
variable s equal next(f)
dump 1 all atom 100 tmp.dump
dump_modify 1 every/time v_s
.. note::
When using a file-style variable with the *every/time* keyword, the
file of timesteps must list a first time that is beyond the time
associated with the current timestep (e.g. it cannot be 0.0). And
it must list one or more times beyond the length of the run you
perform. This is because the dump command will generate an error
if the next time it reads from the file is not a value greater than
the current time. Thus if you wanted output at times 0,15,100 of a
run of length 100 in simulation time, the file should contain the
values 15,100,101 and you should also use the dump_modify first
command. Any final value > 100 could be used in place of 101.
----------
The *first* keyword determines whether a dump snapshot is written on The *first* keyword determines whether a dump snapshot is written on
the very first timestep after the dump command is invoked. This will the very first timestep after the dump command is invoked. This will
always occur if the current timestep is a multiple of N, the frequency always occur if the current timestep is a multiple of $N$, the
specified in the :doc:`dump <dump>` command, including timestep 0. But frequency specified in the :doc:`dump <dump>` command or
if this is not the case, a dump snapshot will only be written if the :doc:`dump_modify every <dump_modify>` command, including timestep 0.
setting of this keyword is *yes*\ . If it is *no*, which is the It will also always occur if the current simulation time is a multiple
of *Delta*, the time interval specified in the doc:`dump_modify
every/time <dump_modify>` command.
But if this is not the case, a dump snapshot will only be written if
the setting of this keyword is *yes*\ . If it is *no*, which is the
default, then it will not be written. default, then it will not be written.
Note that if the argument to the :doc:`dump_modify every
<dump_modify>` doc:`dump_modify every/time <dump_modify>` commands is
a variable and not a numeric value, then specifying *first yes* is the
only way to write a dump snapshot on the first timestep after the dump
command is invoked.
---------- ----------
The *flush* keyword determines whether a flush operation is invoked The *flush* keyword determines whether a flush operation is invoked
@ -342,10 +443,10 @@ The *fileper* keyword is documented below with the *nfile* keyword.
---------- ----------
The *header* keyword toggles whether the dump file will include a header. The *header* keyword toggles whether the dump file will include a
Excluding a header will reduce the size of the dump file for fixes such as header. Excluding a header will reduce the size of the dump file for
:doc:`fix pair/tracker <fix_pair_tracker>` which do not require the information fixes such as :doc:`fix pair/tracker <fix_pair_tracker>` which do not
typically written to the header. require the information typically written to the header.
---------- ----------
@ -641,16 +742,20 @@ threshold criterion is met. Otherwise it is not met.
---------- ----------
The *time* keyword only applies to the dump *atom*, *custom*, and The *time* keyword only applies to the dump *atom*, *custom*, *local*,
*local* styles (and their COMPRESS package versions *atom/gz*, and *xyz* styles (and their COMPRESS package versions *atom/gz*,
*custom/gz* and *local/gz*\ ). If set to *yes*, each frame will will *custom/gz* and *local/gz*\ ). For the first 3 styles, if set to
contain two extra lines before the "ITEM: TIMESTEP" entry: *yes*, each frame will will contain two extra lines before the "ITEM:
TIMESTEP" entry:
.. parsed-literal:: .. parsed-literal::
ITEM: TIME ITEM: TIME
\<elapsed time\> \<elapsed time\>
For the *xyz* style, the simulation time is included on the same line
as the timestep value.
This will output the current elapsed simulation time in current This will output the current elapsed simulation time in current
time units equivalent to the :doc:`thermo keyword <thermo_style>` *time*\ . time units equivalent to the :doc:`thermo keyword <thermo_style>` *time*\ .
This is to simplify post-processing of trajectories using a variable time This is to simplify post-processing of trajectories using a variable time

18
doc/src/fix_dt_reset.rst
View File

@ -78,13 +78,20 @@ outer loop (largest) timestep, which is the same timestep that the
Note that the cumulative simulation time (in time units), which Note that the cumulative simulation time (in time units), which
accounts for changes in the timestep size as a simulation proceeds, accounts for changes in the timestep size as a simulation proceeds,
can be accessed by the :doc:`thermo_style time <thermo_style>` keyword. can be accessed by the :doc:`thermo_style time <thermo_style>`
keyword.
Also note that the :doc:`dump_modify every/time <dump_modify>` option
allows dump files to be written at intervals specified by simulation
time, rather than by timesteps. Simulation time is in time units;
see the :doc:`units <units>` doc page for details.
Restart, fix_modify, output, run start/stop, minimize info Restart, fix_modify, output, run start/stop, minimize info
""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
No information about this fix is written to :doc:`binary restart files <restart>`. None of the :doc:`fix_modify <fix_modify>` options No information about this fix is written to :doc:`binary restart files
are relevant to this fix. <restart>`. None of the :doc:`fix_modify <fix_modify>` options are
relevant to this fix.
This fix computes a global scalar which can be accessed by various This fix computes a global scalar which can be accessed by various
:doc:`output commands <Howto_output>`. The scalar stores the last :doc:`output commands <Howto_output>`. The scalar stores the last
@ -93,7 +100,8 @@ timestep on which the timestep was reset to a new value.
The scalar value calculated by this fix is "intensive". The scalar value calculated by this fix is "intensive".
No parameter of this fix can be used with the *start/stop* keywords of No parameter of this fix can be used with the *start/stop* keywords of
the :doc:`run <run>` command. This fix is not invoked during :doc:`energy minimization <minimize>`. the :doc:`run <run>` command. This fix is not invoked during
:doc:`energy minimization <minimize>`.
Restrictions Restrictions
"""""""""""" """"""""""""
@ -102,7 +110,7 @@ Restrictions
Related commands Related commands
"""""""""""""""" """"""""""""""""
:doc:`timestep <timestep>` :doc:`timestep <timestep>`, :doc:`dump_modify every/time <dump_modify>`
Default Default
""""""" """""""

9
src/COMPRESS/dump_xyz_gz.cpp
View File

@ -88,11 +88,16 @@ void DumpXYZGZ::openfile()
if (multifile) delete[] filecurrent; if (multifile) delete[] filecurrent;
} }
/* ---------------------------------------------------------------------- */
void DumpXYZGZ::write_header(bigint ndump) void DumpXYZGZ::write_header(bigint ndump)
{ {
if (me == 0) { if (me == 0) {
std::string header = fmt::format("{}\n", ndump); auto header = fmt::format("{}\n", ndump);
header += fmt::format("Atoms. Timestep: {}\n", update->ntimestep); if (time_flag) {
double tcurrent = update->atime + (update->ntimestep-update->atimestep) + update->dt;
header += fmt::format(" Atoms. Timestep: {} Time: {:.6f}\n", update->ntimestep, tcurrent);
} else header += fmt::format(" Atoms. Timestep: {}\n", update->ntimestep);
writer.write(header.c_str(), header.length()); writer.write(header.c_str(), header.length());
} }
} }

9
src/COMPRESS/dump_xyz_zstd.cpp
View File

@ -96,11 +96,16 @@ void DumpXYZZstd::openfile()
if (multifile) delete[] filecurrent; if (multifile) delete[] filecurrent;
} }
/* ---------------------------------------------------------------------- */
void DumpXYZZstd::write_header(bigint ndump) void DumpXYZZstd::write_header(bigint ndump)
{ {
if (me == 0) { if (me == 0) {
std::string header = fmt::format("{}\n", ndump); auto header = fmt::format("{}\n", ndump);
header += fmt::format("Atoms. Timestep: {}\n", update->ntimestep); if (time_flag) {
double tcurrent = update->atime + (update->ntimestep-update->atimestep) + update->dt;
header += fmt::format(" Atoms. Timestep: {} Time: {:.6f}\n", update->ntimestep, tcurrent);
} else header += fmt::format(" Atoms. Timestep: {}\n", update->ntimestep);
writer.write(header.c_str(), header.length()); writer.write(header.c_str(), header.length());
} }
} }

8
src/EXTRA-DUMP/dump_dcd.cpp
View File

@ -100,13 +100,17 @@ void DumpDCD::init_style()
if (sort_flag == 0 || sortcol != 0) if (sort_flag == 0 || sortcol != 0)
error->all(FLERR,"Dump dcd requires sorting by atom ID"); error->all(FLERR,"Dump dcd requires sorting by atom ID");
// check that dump frequency has not changed and is not a variable // check that dump modify settings are compatible with dcd
// but only when not being called from the "write_dump" command. // but only when not being called from the "write_dump" command
if (strcmp(id,"WRITE_DUMP") != 0) { if (strcmp(id,"WRITE_DUMP") != 0) {
int idump; int idump;
for (idump = 0; idump < output->ndump; idump++) for (idump = 0; idump < output->ndump; idump++)
if (strcmp(id,output->dump[idump]->id) == 0) break; if (strcmp(id,output->dump[idump]->id) == 0) break;
if (output->mode_dump[idump] == 1)
error->all(FLERR,"Cannot use every/time setting for dump dcd");
if (output->every_dump[idump] == 0) if (output->every_dump[idump] == 0)
error->all(FLERR,"Cannot use variable every setting for dump dcd"); error->all(FLERR,"Cannot use variable every setting for dump dcd");

11
src/EXTRA-DUMP/dump_xtc.cpp
View File

@ -121,17 +121,24 @@ void DumpXTC::init_style()
if (flush_flag) error->all(FLERR,"Cannot set dump_modify flush for dump xtc"); if (flush_flag) error->all(FLERR,"Cannot set dump_modify flush for dump xtc");
// check that dump frequency has not changed and is not a variable // check that dump modify settings are compatible with xtc
// but only when not being called from the "write_dump" command
if (strcmp(id,"WRITE_DUMP") != 0) {
int idump; int idump;
for (idump = 0; idump < output->ndump; idump++) for (idump = 0; idump < output->ndump; idump++)
if (strcmp(id,output->dump[idump]->id) == 0) break; if (strcmp(id,output->dump[idump]->id) == 0) break;
if (output->mode_dump[idump] == 1)
error->all(FLERR,"Cannot use every/time setting for dump xtc");
if (output->every_dump[idump] == 0) if (output->every_dump[idump] == 0)
error->all(FLERR,"Cannot use variable every setting for dump xtc"); error->all(FLERR,"Cannot use every variable setting for dump xtc");
if (nevery_save == 0) nevery_save = output->every_dump[idump]; if (nevery_save == 0) nevery_save = output->every_dump[idump];
else if (nevery_save != output->every_dump[idump]) else if (nevery_save != output->every_dump[idump])
error->all(FLERR,"Cannot change dump_modify every for dump xtc"); error->all(FLERR,"Cannot change dump_modify every for dump xtc");
}
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

30
src/dump.cpp
View File

@ -918,11 +918,6 @@ void Dump::modify_params(int narg, char **arg)
else delay_flag = 0; else delay_flag = 0;
iarg += 2; iarg += 2;
} else if (strcmp(arg[iarg],"header") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
header_flag = utils::logical(FLERR,arg[iarg+1],false,lmp);
iarg += 2;
} else if (strcmp(arg[iarg],"every") == 0) { } else if (strcmp(arg[iarg],"every") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
int idump; int idump;
@ -937,9 +932,28 @@ void Dump::modify_params(int narg, char **arg)
n = utils::inumeric(FLERR,arg[iarg+1],false,lmp); n = utils::inumeric(FLERR,arg[iarg+1],false,lmp);
if (n <= 0) error->all(FLERR,"Illegal dump_modify command"); if (n <= 0) error->all(FLERR,"Illegal dump_modify command");
} }
output->mode_dump[idump] = 0;
output->every_dump[idump] = n; output->every_dump[idump] = n;
iarg += 2; iarg += 2;
} else if (strcmp(arg[iarg],"every/time") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
int idump;
for (idump = 0; idump < output->ndump; idump++)
if (strcmp(id,output->dump[idump]->id) == 0) break;
double delta;
if (strstr(arg[iarg+1],"v_") == arg[iarg+1]) {
delete [] output->var_dump[idump];
output->var_dump[idump] = utils::strdup(&arg[iarg+1][2]);
delta = 0.0;
} else {
delta = utils::numeric(FLERR,arg[iarg+1],false,lmp);
if (delta <= 0.0) error->all(FLERR,"Illegal dump_modify command");
}
output->mode_dump[idump] = 1;
output->every_time_dump[idump] = delta;
iarg += 2;
} else if (strcmp(arg[iarg],"fileper") == 0) { } else if (strcmp(arg[iarg],"fileper") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
if (!multiproc) if (!multiproc)
@ -1007,6 +1021,11 @@ void Dump::modify_params(int narg, char **arg)
iarg += n; iarg += n;
} }
} else if (strcmp(arg[iarg],"header") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
header_flag = utils::logical(FLERR,arg[iarg+1],false,lmp);
iarg += 2;
} else if (strcmp(arg[iarg],"maxfiles") == 0) { } else if (strcmp(arg[iarg],"maxfiles") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command");
if (!multifile) if (!multifile)
@ -1113,6 +1132,7 @@ double Dump::compute_time()
{ {
return update->atime + (update->ntimestep - update->atimestep)*update->dt; return update->atime + (update->ntimestep - update->atimestep)*update->dt;
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
return # of bytes of allocated memory return # of bytes of allocated memory
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */

4
src/dump.h
View File

@ -26,7 +26,7 @@ class Dump : protected Pointers {
int igroup, groupbit; // group that Dump is performed on int igroup, groupbit; // group that Dump is performed on
int first_flag; // 0 if no initial dump, 1 if yes initial dump int first_flag; // 0 if no initial dump, 1 if yes initial dump
int clearstep; // 1 if dump invokes computes, 0 if not int clearstep; // 1 if dump can invoke computes, 0 if not
int comm_forward; // size of forward communication (0 if none) int comm_forward; // size of forward communication (0 if none)
int comm_reverse; // size of reverse communication (0 if none) int comm_reverse; // size of reverse communication (0 if none)
@ -75,7 +75,7 @@ class Dump : protected Pointers {
int sortcol; // 0 to sort on ID, 1-N on columns int sortcol; // 0 to sort on ID, 1-N on columns
int sortcolm1; // sortcol - 1 int sortcolm1; // sortcol - 1
int sortorder; // ASCEND or DESCEND int sortorder; // ASCEND or DESCEND
int time_flag; // 1 if output accumulated time int time_flag; // 1 if output simulation time
int unit_flag; // 1 if dump should contain unit information int unit_flag; // 1 if dump should contain unit information
int unit_count; // # of times the unit information was written int unit_count; // # of times the unit information was written
int delay_flag; // 1 if delay output until delaystep int delay_flag; // 1 if delay output until delaystep

7
src/dump_xyz.cpp
View File

@ -130,8 +130,10 @@ int DumpXYZ::modify_param(int narg, char **arg)
void DumpXYZ::write_header(bigint n) void DumpXYZ::write_header(bigint n)
{ {
if (me == 0) { if (me == 0) {
fprintf(fp,BIGINT_FORMAT "\n",n); if (time_flag) {
fprintf(fp,"Atoms. Timestep: " BIGINT_FORMAT "\n",update->ntimestep); double tcurrent = update->atime + (update->ntimestep-update->atimestep) + update->dt;
fmt::print(fp,"{}\n Atoms. Timestep: {} Time: {:.6f}\n", n, update->ntimestep, tcurrent);
} else fmt::print(fp,"{}\n Atoms. Timestep: {}\n", n, update->ntimestep);
} }
} }
@ -159,7 +161,6 @@ void DumpXYZ::pack(tagint *ids)
} }
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
convert mybuf of doubles to one big formatted string in sbuf convert mybuf of doubles to one big formatted string in sbuf
return -1 if strlen exceeds an int, since used as arg in MPI calls in Dump return -1 if strlen exceeds an int, since used as arg in MPI calls in Dump

12
src/fix_dt_reset.cpp
View File

@ -121,14 +121,6 @@ void FixDtReset::init()
respaflag = 0; respaflag = 0;
if (utils::strmatch(update->integrate_style, "^respa")) respaflag = 1; if (utils::strmatch(update->integrate_style, "^respa")) respaflag = 1;
// check for DCD or XTC dumps
for (int i = 0; i < output->ndump; i++)
if ((strcmp(output->dump[i]->style, "dcd") == 0 ||
strcmp(output->dump[i]->style, "xtc") == 0) &&
comm->me == 0)
error->warning(FLERR, "Dump dcd/xtc timestamp may be wrong with fix dt/reset");
ftm2v = force->ftm2v; ftm2v = force->ftm2v;
mvv2e = force->mvv2e; mvv2e = force->mvv2e;
dt = update->dt; dt = update->dt;
@ -197,12 +189,16 @@ void FixDtReset::end_of_step()
laststep = update->ntimestep; laststep = update->ntimestep;
// calls to other classes that need to know timestep size changed
// similar logic is in Input::timestep()
update->update_time(); update->update_time();
update->dt = dt; update->dt = dt;
update->dt_default = 0; update->dt_default = 0;
if (respaflag) update->integrate->reset_dt(); if (respaflag) update->integrate->reset_dt();
if (force->pair) force->pair->reset_dt(); if (force->pair) force->pair->reset_dt();
for (int i = 0; i < modify->nfix; i++) modify->fix[i]->reset_dt(); for (int i = 0; i < modify->nfix; i++) modify->fix[i]->reset_dt();
output->reset_dt();
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

19
src/input.cpp
View File

@ -27,9 +27,11 @@
#include "dihedral.h" #include "dihedral.h"
#include "domain.h" #include "domain.h"
#include "error.h" #include "error.h"
#include "fix.h"
#include "force.h" #include "force.h"
#include "group.h" #include "group.h"
#include "improper.h" #include "improper.h"
#include "integrate.h"
#include "kspace.h" #include "kspace.h"
#include "memory.h" #include "memory.h"
#include "min.h" #include "min.h"
@ -1836,8 +1838,25 @@ void Input::timer_command()
void Input::timestep() void Input::timestep()
{ {
if (narg != 1) error->all(FLERR,"Illegal timestep command"); if (narg != 1) error->all(FLERR,"Illegal timestep command");
update->update_time();
update->dt = utils::numeric(FLERR,arg[0],false,lmp); update->dt = utils::numeric(FLERR,arg[0],false,lmp);
update->dt_default = 0; update->dt_default = 0;
// timestep command can be invoked between runs or by run every
// calls to other classes that need to know timestep size changed
// similar logic is in FixDtReset::end_of_step()
// only need to do this if a run has already occurred
if (update->first_update == 0) return;
int respaflag = 0;
if (utils::strmatch(update->integrate_style, "^respa")) respaflag = 1;
if (respaflag) update->integrate->reset_dt();
if (force->pair) force->pair->reset_dt();
for (int i = 0; i < modify->nfix; i++) modify->fix[i]->reset_dt();
output->reset_dt();
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

18
src/integrate.cpp
View File

@ -20,6 +20,7 @@
#include "kspace.h" #include "kspace.h"
#include "modify.h" #include "modify.h"
#include "pair.h" #include "pair.h"
#include "output.h"
#include "update.h" #include "update.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
@ -115,7 +116,12 @@ void Integrate::ev_setup()
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
set eflag,vflag for current iteration set eflag,vflag for current iteration
based on computes that need energy/virial info on this timestep based on
(1) computes that need energy/virial info on this timestep
(2) time dumps that may need per-atom compute info on this timestep
NOTE: inefficient to add all per-atom eng/virial computes
but don't know which ones the dump needs
see NOTE in output.cpp
invoke matchstep() on all timestep-dependent computes to clear their arrays invoke matchstep() on all timestep-dependent computes to clear their arrays
eflag: set any or no bits eflag: set any or no bits
ENERGY_GLOBAL bit for global energy ENERGY_GLOBAL bit for global energy
@ -133,6 +139,10 @@ void Integrate::ev_set(bigint ntimestep)
{ {
int i,flag; int i,flag;
int tdflag = 0;
if (output->any_time_dumps &&
output->next_time_dump_any == ntimestep) tdflag = 1;
flag = 0; flag = 0;
int eflag_global = 0; int eflag_global = 0;
for (i = 0; i < nelist_global; i++) for (i = 0; i < nelist_global; i++)
@ -143,7 +153,7 @@ void Integrate::ev_set(bigint ntimestep)
int eflag_atom = 0; int eflag_atom = 0;
for (i = 0; i < nelist_atom; i++) for (i = 0; i < nelist_atom; i++)
if (elist_atom[i]->matchstep(ntimestep)) flag = 1; if (elist_atom[i]->matchstep(ntimestep)) flag = 1;
if (flag) eflag_atom = ENERGY_ATOM; if (flag || (tdflag && nelist_atom)) eflag_atom = ENERGY_ATOM;
if (eflag_global) update->eflag_global = ntimestep; if (eflag_global) update->eflag_global = ntimestep;
if (eflag_atom) update->eflag_atom = ntimestep; if (eflag_atom) update->eflag_atom = ntimestep;
@ -159,13 +169,13 @@ void Integrate::ev_set(bigint ntimestep)
int vflag_atom = 0; int vflag_atom = 0;
for (i = 0; i < nvlist_atom; i++) for (i = 0; i < nvlist_atom; i++)
if (vlist_atom[i]->matchstep(ntimestep)) flag = 1; if (vlist_atom[i]->matchstep(ntimestep)) flag = 1;
if (flag) vflag_atom = VIRIAL_ATOM; if (flag || (tdflag && nvlist_atom)) vflag_atom = VIRIAL_ATOM;
flag = 0; flag = 0;
int cvflag_atom = 0; int cvflag_atom = 0;
for (i = 0; i < ncvlist_atom; i++) for (i = 0; i < ncvlist_atom; i++)
if (cvlist_atom[i]->matchstep(ntimestep)) flag = 1; if (cvlist_atom[i]->matchstep(ntimestep)) flag = 1;
if (flag) cvflag_atom = VIRIAL_CENTROID; if (flag || (tdflag && ncvlist_atom)) cvflag_atom = VIRIAL_CENTROID;
if (vflag_global) update->vflag_global = ntimestep; if (vflag_global) update->vflag_global = ntimestep;
if (vflag_atom || cvflag_atom) update->vflag_atom = ntimestep; if (vflag_atom || cvflag_atom) update->vflag_atom = ntimestep;

399
src/output.cpp
View File

@ -12,6 +12,10 @@
See the README file in the top-level LAMMPS directory. See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Michal Kanski (Jagiellonian U) for simulation time dumps
------------------------------------------------------------------------- */
#include "output.h" #include "output.h"
#include "style_dump.h" // IWYU pragma: keep #include "style_dump.h" // IWYU pragma: keep
@ -29,11 +33,13 @@
#include "variable.h" #include "variable.h"
#include "write_restart.h" #include "write_restart.h"
#include <cmath>
#include <cstring> #include <cstring>
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
#define DELTA 1 #define DELTA 1
#define EPSDT 1.0e-6
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
initialize all output initialize all output
@ -59,8 +65,11 @@ Output::Output(LAMMPS *lmp) : Pointers(lmp)
ndump = 0; ndump = 0;
max_dump = 0; max_dump = 0;
mode_dump = nullptr;
every_dump = nullptr; every_dump = nullptr;
every_time_dump = nullptr;
next_dump = nullptr; next_dump = nullptr;
next_time_dump = nullptr;
last_dump = nullptr; last_dump = nullptr;
var_dump = nullptr; var_dump = nullptr;
ivar_dump = nullptr; ivar_dump = nullptr;
@ -92,8 +101,11 @@ Output::~Output()
if (thermo) delete thermo; if (thermo) delete thermo;
delete [] var_thermo; delete [] var_thermo;
memory->destroy(mode_dump);
memory->destroy(every_dump); memory->destroy(every_dump);
memory->destroy(every_time_dump);
memory->destroy(next_dump); memory->destroy(next_dump);
memory->destroy(next_time_dump);
memory->destroy(last_dump); memory->destroy(last_dump);
for (int i = 0; i < ndump; i++) delete [] var_dump[i]; for (int i = 0; i < ndump; i++) delete [] var_dump[i];
memory->sfree(var_dump); memory->sfree(var_dump);
@ -125,13 +137,17 @@ void Output::init()
} }
for (int i = 0; i < ndump; i++) dump[i]->init(); for (int i = 0; i < ndump; i++) dump[i]->init();
for (int i = 0; i < ndump; i++) any_time_dumps = 0;
if (every_dump[i] == 0) { for (int i = 0; i < ndump; i++) {
if (mode_dump[i]) any_time_dumps = 1;
if ((mode_dump[i] == 0 && every_dump[i] == 0) ||
(mode_dump[i] == 1 && every_time_dump[i] == 0.0)) {
ivar_dump[i] = input->variable->find(var_dump[i]); ivar_dump[i] = input->variable->find(var_dump[i]);
if (ivar_dump[i] < 0) if (ivar_dump[i] < 0)
error->all(FLERR,"Variable name for dump every does not exist"); error->all(FLERR,"Variable name for dump every or delta does not exist");
if (!input->variable->equalstyle(ivar_dump[i])) if (!input->variable->equalstyle(ivar_dump[i]))
error->all(FLERR,"Variable for dump every is invalid style"); error->all(FLERR,"Variable for dump every or delta is invalid style");
}
} }
if (restart_flag_single && restart_every_single == 0) { if (restart_flag_single && restart_every_single == 0) {
@ -161,49 +177,80 @@ void Output::setup(int memflag)
{ {
bigint ntimestep = update->ntimestep; bigint ntimestep = update->ntimestep;
// perform dump at start of run only if: // consider all dumps
// current timestep is multiple of every and last dump not >= this step // decide whether to write snapshot and/or calculate next step for dump
// this is first run after dump created and firstflag is set
// note that variable freq will not write unless triggered by firstflag
// set next_dump to multiple of every or variable value
// set next_dump_any to smallest next_dump
// wrap dumps that invoke computes and variable eval with clear/add
// if dump not written now, use addstep_compute_all() since don't know
// what computes the dump write would invoke
// if no dumps, set next_dump_any to last+1 so will not influence next
int writeflag;
if (ndump && update->restrict_output == 0) { if (ndump && update->restrict_output == 0) {
next_time_dump_any = MAXBIGINT;
for (int idump = 0; idump < ndump; idump++) { for (int idump = 0; idump < ndump; idump++) {
if (dump[idump]->clearstep || every_dump[idump] == 0)
// wrap step dumps that invoke computes or do variable eval with clear/add
// see NOTE in write() about also wrapping time dumps
if (mode_dump[idump] == 0 &&
(dump[idump]->clearstep || var_dump[idump]))
modify->clearstep_compute(); modify->clearstep_compute();
writeflag = 0;
if (every_dump[idump] && ntimestep % every_dump[idump] == 0 && // write a snapshot at setup only if any of these 3 conditions hold
last_dump[idump] != ntimestep) writeflag = 1; // (1) this is first run since dump was created and its first_flag = 0
// (2) mode_dump = 0 and timestep is multiple of every_dump
// (3) mode_dump = 1 and time is multiple of every_time_dump (within EPSDT)
// (2) and (3) only apply for non-variable dump intervals
// finally, do not write if same snapshot written previously,
// i.e. on last timestep of previous run
int writeflag = 0;
if (last_dump[idump] < 0 && dump[idump]->first_flag == 1) writeflag = 1; if (last_dump[idump] < 0 && dump[idump]->first_flag == 1) writeflag = 1;
if (mode_dump[idump] == 0) {
if (every_dump[idump] && (ntimestep % every_dump[idump] == 0))
writeflag = 1;
} else {
if (every_time_dump[idump] > 0.0) {
double tcurrent = update->atime +
(ntimestep - update->atimestep) * update->dt;
double remainder = fmod(tcurrent,every_time_dump[idump]);
if ((remainder < EPSDT*update->dt) ||
(every_time_dump[idump] - remainder < EPSDT*update->dt))
writeflag = 1;
}
}
if (last_dump[idump] == ntimestep) writeflag = 0;
// perform dump
if (writeflag) { if (writeflag) {
dump[idump]->write(); dump[idump]->write();
last_dump[idump] = ntimestep; last_dump[idump] = ntimestep;
} }
if (every_dump[idump])
next_dump[idump] = // calculate timestep and/or time for next dump
(ntimestep/every_dump[idump])*every_dump[idump] + every_dump[idump]; // set next_dump and next_time_dump, 0 arg for setup()
else { // only do this if dump written or dump has not been written yet
bigint nextdump = static_cast<bigint>
(input->variable->compute_equal(ivar_dump[idump])); if (writeflag || last_dump[idump] < 0)
if (nextdump <= ntimestep) calculate_next_dump(0,idump,ntimestep);
error->all(FLERR,"Dump every variable returned a bad timestep");
next_dump[idump] = nextdump; // if dump not written now, use addstep_compute_all()
} // since don't know what computes the dump will invoke
if (dump[idump]->clearstep || every_dump[idump] == 0) {
if (mode_dump[idump] == 0 &&
(dump[idump]->clearstep || var_dump[idump])) {
if (writeflag) modify->addstep_compute(next_dump[idump]); if (writeflag) modify->addstep_compute(next_dump[idump]);
else modify->addstep_compute_all(next_dump[idump]); else modify->addstep_compute_all(next_dump[idump]);
} }
if (mode_dump[idump] && (dump[idump]->clearstep || var_dump[idump]))
next_time_dump_any = MIN(next_time_dump_any,next_dump[idump]);
if (idump) next_dump_any = MIN(next_dump_any,next_dump[idump]); if (idump) next_dump_any = MIN(next_dump_any,next_dump[idump]);
else next_dump_any = next_dump[0]; else next_dump_any = next_dump[0];
} }
// if no dumps, set next_dump_any to last+1 so will not influence next
} else next_dump_any = update->laststep + 1; } else next_dump_any = update->laststep + 1;
// do not write restart files at start of run // do not write restart files at start of run
@ -271,39 +318,60 @@ void Output::setup(int memflag)
next = MIN(next_dump_any,next_restart); next = MIN(next_dump_any,next_restart);
next = MIN(next,next_thermo); next = MIN(next,next_thermo);
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
perform all output for this timestep perform all output for this timestep
only perform output if next matches current step and last output doesn't only perform output if next matches current step and last output doesn't
do dump/restart before thermo so thermo CPU time will include them do dump/restart before thermo so thermo CPU time will include them
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Output::write(bigint ntimestep) void Output::write(bigint ntimestep)
{ {
// next_dump does not force output on last step of run // perform dump if its next_dump = current ntimestep
// wrap dumps that invoke computes or eval of variable with clear/add // but not if it was already written on this step
// set next_dump and also next_time_dump for mode_dump = 1
// set next_dump_any to smallest next_dump
// wrap step dumps that invoke computes or do variable eval with clear/add
// NOTE:
// not wrapping time dumps means that Integrate::ev_set()
// needs to trigger all per-atom eng/virial computes
// on a timestep where any time dump will be output
// could wrap time dumps as well, if timestep size did not vary
// if wrap when timestep size varies frequently,
// then can do many unneeded addstep() --> inefficient
// hard to know if timestep varies, since run every could change it
// can't remove an uneeded addstep from a compute, b/c don't know
// what other command may have added it
int writeflag;
if (next_dump_any == ntimestep) { if (next_dump_any == ntimestep) {
for (int idump = 0; idump < ndump; idump++) { for (int idump = 0; idump < ndump; idump++) {
next_time_dump_any = MAXBIGINT;
if (next_dump[idump] == ntimestep) { if (next_dump[idump] == ntimestep) {
if (dump[idump]->clearstep || every_dump[idump] == 0) if (last_dump[idump] == ntimestep) continue;
if (mode_dump[idump] == 0 &&
(dump[idump]->clearstep || var_dump[idump]))
modify->clearstep_compute(); modify->clearstep_compute();
if (last_dump[idump] != ntimestep) {
// perform dump
// reset next_dump and next_time_dump, 1 arg for write()
dump[idump]->write(); dump[idump]->write();
last_dump[idump] = ntimestep; last_dump[idump] = ntimestep;
} calculate_next_dump(1,idump,ntimestep);
if (every_dump[idump]) next_dump[idump] += every_dump[idump];
else { if (mode_dump[idump] == 0 &&
bigint nextdump = static_cast<bigint> (dump[idump]->clearstep || var_dump[idump]))
(input->variable->compute_equal(ivar_dump[idump]));
if (nextdump <= ntimestep)
error->all(FLERR,"Dump every variable returned a bad timestep");
next_dump[idump] = nextdump;
}
if (dump[idump]->clearstep || every_dump[idump] == 0)
modify->addstep_compute(next_dump[idump]); modify->addstep_compute(next_dump[idump]);
} }
if (mode_dump[idump] && (dump[idump]->clearstep || var_dump[idump]))
next_time_dump_any = MIN(next_time_dump_any,next_dump[idump]);
if (idump) next_dump_any = MIN(next_dump_any,next_dump[idump]); if (idump) next_dump_any = MIN(next_dump_any,next_dump[idump]);
else next_dump_any = next_dump[0]; else next_dump_any = next_dump[0];
} }
@ -311,7 +379,7 @@ void Output::write(bigint ntimestep)
// next_restart does not force output on last step of run // next_restart does not force output on last step of run
// for toggle = 0, replace "*" with current timestep in restart filename // for toggle = 0, replace "*" with current timestep in restart filename
// eval of variable may invoke computes so wrap with clear/add // next restart variable may invoke computes so wrap with clear/add
if (next_restart == ntimestep) { if (next_restart == ntimestep) {
if (next_restart_single == ntimestep) { if (next_restart_single == ntimestep) {
@ -381,28 +449,141 @@ void Output::write(bigint ntimestep)
next = MIN(next_dump_any,next_restart); next = MIN(next_dump_any,next_restart);
next = MIN(next,next_thermo); next = MIN(next,next_thermo);
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
force a snapshot to be written for all dumps force a snapshot to be written for all dumps
called from PRD and TAD called from PRD and TAD
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Output::write_dump(bigint ntimestep) void Output::write_dump(bigint ntimestep)
{ {
for (int idump = 0; idump < ndump; idump++) { for (int idump = 0; idump < ndump; idump++) {
dump[idump]->write(); dump[idump]->write();
last_dump[idump] = ntimestep; last_dump[idump] = ntimestep;
} }
}
/* ----------------------------------------------------------------------
calculate when next dump occurs for Dump instance idump
operates in one of two modes, based on mode_dump flag
for timestep mode, set next_dump
for simulation time mode, set next_time_dump and next_dump
which flag depends on caller
0 = from setup() at start of run
1 = from write() during run each time a dump file is written
2 = from reset_dt() called from fix dt/reset when it changes timestep size
------------------------------------------------------------------------- */
void Output::calculate_next_dump(int which, int idump, bigint ntimestep)
{
// dump mode is by timestep
// just set next_dump
if (mode_dump[idump] == 0) {
if (every_dump[idump]) {
// which = 0: nextdump = next multiple of every_dump
// which = 1: increment nextdump by every_dump
if (which == 0)
next_dump[idump] =
(ntimestep/every_dump[idump])*every_dump[idump] + every_dump[idump];
else if (which == 1)
next_dump[idump] += every_dump[idump];
} else {
next_dump[idump] = static_cast<bigint>
(input->variable->compute_equal(ivar_dump[idump]));
if (next_dump[idump] <= ntimestep)
error->all(FLERR,"Dump every variable returned a bad timestep");
}
// dump mode is by simulation time
// set next_time_dump and next_dump
} else {
bigint nextdump;
double nexttime;
double tcurrent = update->atime +
(ntimestep - update->atimestep) * update->dt;
if (every_time_dump[idump] > 0.0) {
// which = 0: nexttime = next multiple of every_time_dump
// which = 1: increment nexttime by every_time_dump
// which = 2: no change to previous nexttime (only timestep has changed)
if (which == 0)
nexttime = static_cast<bigint> (tcurrent/every_time_dump[idump]) *
every_time_dump[idump] + every_time_dump[idump];
else if (which == 1)
nexttime = next_time_dump[idump] + every_time_dump[idump];
else if (which == 2)
nexttime = next_time_dump[idump];
nextdump = ntimestep +
static_cast<bigint> ((nexttime - tcurrent - EPSDT*update->dt) /
update->dt) + 1;
// if delta is too small to reach next timestep, use multiple of delta
if (nextdump == ntimestep) {
double tnext = update->atime +
(ntimestep+1 - update->atimestep) * update->dt;
int multiple = static_cast<int>
((tnext - nexttime) / every_time_dump[idump]);
nexttime = nexttime + (multiple+1)*every_time_dump[idump];
nextdump = ntimestep +
static_cast<bigint> ((nexttime - tcurrent - EPSDT*update->dt) /
update->dt) + 1;
}
} else {
// do not re-evaulate variable for which = 2, leave nexttime as-is
// unless next_time_dump < 0.0, which means variable never yet evaluated
if (which < 2 || next_time_dump[idump] < 0.0) {
nexttime = input->variable->compute_equal(ivar_dump[idump]);
} else
nexttime = next_time_dump[idump];
if (nexttime <= tcurrent)
error->all(FLERR,"Dump every/time variable returned a bad time");
nextdump = ntimestep +
static_cast<bigint> ((nexttime - tcurrent - EPSDT*update->dt) /
update->dt) + 1;
if (nextdump <= ntimestep)
error->all(FLERR,"Dump every/time variable too small for next timestep");
}
next_time_dump[idump] = nexttime;
next_dump[idump] = nextdump;
}
}
/* ---------------------------------------------------------------------- */
int Output::check_time_dumps(bigint ntimestep)
{
int nowflag = 0;
for (int i = 0; i < ndump; i++)
if (mode_dump[i] && next_dump[i] == ntimestep) nowflag = 1;
return nowflag;
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
force restart file(s) to be written force restart file(s) to be written
called from PRD and TAD called from PRD and TAD
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Output::write_restart(bigint ntimestep) void Output::write_restart(bigint ntimestep)
{ {
if (restart_flag_single) { if (restart_flag_single) {
std::string file = restart1; std::string file = restart1;
std::size_t found = file.find('*'); std::size_t found = file.find('*');
@ -422,46 +603,27 @@ void Output::write_restart(bigint ntimestep)
} }
last_restart = ntimestep; last_restart = ntimestep;
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
timestep is being changed, called by update->reset_timestep() timestep is being changed, called by update->reset_timestep()
reset next timestep values for dumps, restart, thermo output for dumps, require that no dump is "active"
meaning that a snapshot has already been output
reset next output values for restart and thermo
reset to smallest value >= new timestep reset to smallest value >= new timestep
if next timestep set by variable evaluation, if next timestep set by variable evaluation,
eval for ntimestep-1, so current ntimestep can be returned if needed eval for ntimestep-1, so current ntimestep can be returned if needed
no guarantee that variable can be evaluated for ntimestep-1 no guarantee that variable can be evaluated for ntimestep-1
if it depends on computes, but live with that rare case for now e.g. if it depends on computes, but live with that rare case for now
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Output::reset_timestep(bigint ntimestep) void Output::reset_timestep(bigint ntimestep)
{ {
next_dump_any = MAXBIGINT; next_dump_any = MAXBIGINT;
for (int idump = 0; idump < ndump; idump++) { for (int idump = 0; idump < ndump; idump++)
if (every_dump[idump]) { if (last_dump[idump] >= 0)
next_dump[idump] = (ntimestep/every_dump[idump])*every_dump[idump]; error->all(FLERR,
if (next_dump[idump] < ntimestep) next_dump[idump] += every_dump[idump]; "Cannot reset timestep with active dump - must undump first");
} else {
// ivar_dump may not be initialized
if (ivar_dump[idump] < 0) {
ivar_dump[idump] = input->variable->find(var_dump[idump]);
if (ivar_dump[idump] < 0)
error->all(FLERR,"Variable name for dump every does not exist");
if (!input->variable->equalstyle(ivar_dump[idump]))
error->all(FLERR,"Variable for dump every is invalid style");
}
modify->clearstep_compute();
update->ntimestep--;
bigint nextdump = static_cast<bigint>
(input->variable->compute_equal(ivar_dump[idump]));
if (nextdump < ntimestep)
error->all(FLERR,"Dump every variable returned a bad timestep");
update->ntimestep++;
next_dump[idump] = nextdump;
modify->addstep_compute(next_dump[idump]);
}
next_dump_any = MIN(next_dump_any,next_dump[idump]);
}
if (restart_flag_single) { if (restart_flag_single) {
if (restart_every_single) { if (restart_every_single) {
@ -519,10 +681,43 @@ void Output::reset_timestep(bigint ntimestep)
next_thermo = MIN(next_thermo,update->laststep); next_thermo = MIN(next_thermo,update->laststep);
} else next_thermo = update->laststep; } else next_thermo = update->laststep;
next = MIN(next_dump_any,next_restart);
next = MIN(next,next_thermo);
}
/* ----------------------------------------------------------------------
timestep size is being changed
reset next output values for dumps which have mode_dump=1
called by fix dt/reset (at end of step)
or called by timestep command via run every (also at end of step)
------------------------------------------------------------------------- */
void Output::reset_dt()
{
bigint ntimestep = update->ntimestep;
next_time_dump_any = MAXBIGINT;
for (int idump = 0; idump < ndump; idump++) {
if (mode_dump[idump] == 0) continue;
// reset next_dump but do not change next_time_dump, 2 arg for reset_dt()
// do not invoke for a dump already scheduled for this step
// since timestep change affects next step
if (next_dump[idump] != ntimestep)
calculate_next_dump(2,idump,update->ntimestep);
if (dump[idump]->clearstep || var_dump[idump])
next_time_dump_any = MIN(next_time_dump_any,next_dump[idump]);
}
next_dump_any = MIN(next_dump_any,next_time_dump_any);
next = MIN(next_dump_any,next_restart); next = MIN(next_dump_any,next_restart);
next = MIN(next,next_thermo); next = MIN(next,next_thermo);
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
add a Dump to list of Dumps add a Dump to list of Dumps
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
@ -547,21 +742,17 @@ void Output::add_dump(int narg, char **arg)
max_dump += DELTA; max_dump += DELTA;
dump = (Dump **) dump = (Dump **)
memory->srealloc(dump,max_dump*sizeof(Dump *),"output:dump"); memory->srealloc(dump,max_dump*sizeof(Dump *),"output:dump");
memory->grow(mode_dump,max_dump,"output:mode_dump");
memory->grow(every_dump,max_dump,"output:every_dump"); memory->grow(every_dump,max_dump,"output:every_dump");
memory->grow(every_time_dump,max_dump,"output:every_time_dump");
memory->grow(next_dump,max_dump,"output:next_dump"); memory->grow(next_dump,max_dump,"output:next_dump");
memory->grow(next_time_dump,max_dump,"output:next_time_dump");
memory->grow(last_dump,max_dump,"output:last_dump"); memory->grow(last_dump,max_dump,"output:last_dump");
var_dump = (char **) var_dump = (char **)
memory->srealloc(var_dump,max_dump*sizeof(char *),"output:var_dump"); memory->srealloc(var_dump,max_dump*sizeof(char *),"output:var_dump");
memory->grow(ivar_dump,max_dump,"output:ivar_dump"); memory->grow(ivar_dump,max_dump,"output:ivar_dump");
} }
// initialize per-dump data to suitable default values
every_dump[ndump] = 0;
last_dump[ndump] = -1;
var_dump[ndump] = nullptr;
ivar_dump[ndump] = -1;
// create the Dump // create the Dump
if (dump_map->find(arg[2]) != dump_map->end()) { if (dump_map->find(arg[2]) != dump_map->end()) {
@ -569,10 +760,17 @@ void Output::add_dump(int narg, char **arg)
dump[ndump] = dump_creator(lmp, narg, arg); dump[ndump] = dump_creator(lmp, narg, arg);
} else error->all(FLERR,utils::check_packages_for_style("dump",arg[2],lmp)); } else error->all(FLERR,utils::check_packages_for_style("dump",arg[2],lmp));
// initialize per-dump data to suitable default values
mode_dump[ndump] = 0;
every_dump[ndump] = utils::inumeric(FLERR,arg[3],false,lmp); every_dump[ndump] = utils::inumeric(FLERR,arg[3],false,lmp);
if (every_dump[ndump] <= 0) error->all(FLERR,"Illegal dump command"); if (every_dump[ndump] <= 0) error->all(FLERR,"Illegal dump command");
every_time_dump[ndump] = 0.0;
next_time_dump[ndump] = -1.0;
last_dump[ndump] = -1; last_dump[ndump] = -1;
var_dump[ndump] = nullptr; var_dump[ndump] = nullptr;
ivar_dump[ndump] = -1;
ndump++; ndump++;
} }
@ -624,8 +822,11 @@ void Output::delete_dump(char *id)
for (int i = idump+1; i < ndump; i++) { for (int i = idump+1; i < ndump; i++) {
dump[i-1] = dump[i]; dump[i-1] = dump[i];
mode_dump[i-1] = mode_dump[i];
every_dump[i-1] = every_dump[i]; every_dump[i-1] = every_dump[i];
every_time_dump[i-1] = every_time_dump[i];
next_dump[i-1] = next_dump[i]; next_dump[i-1] = next_dump[i];
next_time_dump[i-1] = next_time_dump[i];
last_dump[i-1] = last_dump[i]; last_dump[i-1] = last_dump[i];
var_dump[i-1] = var_dump[i]; var_dump[i-1] = var_dump[i];
ivar_dump[i-1] = ivar_dump[i]; ivar_dump[i-1] = ivar_dump[i];

20
src/output.h
View File

@ -35,12 +35,17 @@ class Output : protected Pointers {
int ndump; // # of Dumps defined int ndump; // # of Dumps defined
int max_dump; // max size of Dump list int max_dump; // max size of Dump list
bigint next_dump_any; // next timestep for any Dump bigint next_dump_any; // next timestep for any dump
int *every_dump; // write freq for each Dump, 0 if var bigint next_time_dump_any; // next timestep for any time dump with computes
bigint *next_dump; // next timestep to do each Dump int any_time_dumps; // 1 if any time dump defined
int *mode_dump; // 0/1 if write every N timesteps or Delta in sim time
int *every_dump; // dump every N timesteps, 0 if variable
double *every_time_dump; // dump every Delta of sim time, 0.0 if variable
bigint *next_dump; // next timestep to perform dump
double *next_time_dump; // next simulation time to perform dump (mode = 1)
bigint *last_dump; // last timestep each snapshot was output bigint *last_dump; // last timestep each snapshot was output
char **var_dump; // variable name for dump frequency char **var_dump; // variable name for next dump (steps or sim time)
int *ivar_dump; // variable index for dump frequency int *ivar_dump; // variable index of var_dump name
Dump **dump; // list of defined Dumps Dump **dump; // list of defined Dumps
int restart_flag; // 1 if any restart files are written int restart_flag; // 1 if any restart files are written
@ -72,12 +77,14 @@ class Output : protected Pointers {
void write(bigint); // output for current timestep void write(bigint); // output for current timestep
void write_dump(bigint); // force output of dump snapshots void write_dump(bigint); // force output of dump snapshots
void write_restart(bigint); // force output of a restart file void write_restart(bigint); // force output of a restart file
void reset_timestep(bigint); // reset next timestep for all output void reset_timestep(bigint); // reset output which depends on timestep
void reset_dt(); // reset output which depends on timestep size
void add_dump(int, char **); // add a Dump to Dump list void add_dump(int, char **); // add a Dump to Dump list
void modify_dump(int, char **); // modify a Dump void modify_dump(int, char **); // modify a Dump
void delete_dump(char *); // delete a Dump from Dump list void delete_dump(char *); // delete a Dump from Dump list
int find_dump(const char *); // find a Dump ID int find_dump(const char *); // find a Dump ID
int check_time_dumps(bigint); // check if any time dump is output now
void set_thermo(int, char **); // set thermo output freqquency void set_thermo(int, char **); // set thermo output freqquency
void create_thermo(int, char **); // create a thermo style void create_thermo(int, char **); // create a thermo style
@ -87,6 +94,7 @@ class Output : protected Pointers {
private: private:
template <typename T> static Dump *dump_creator(LAMMPS *, int, char **); template <typename T> static Dump *dump_creator(LAMMPS *, int, char **);
void calculate_next_dump(int, int, bigint);
}; };
} // namespace LAMMPS_NS } // namespace LAMMPS_NS

27
src/update.cpp
View File

@ -458,7 +458,7 @@ Min *Update::minimize_creator(LAMMPS *lmp)
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
reset timestep as called from input script reset timestep called from input script
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Update::reset_timestep(int narg, char **arg) void Update::reset_timestep(int narg, char **arg)
@ -470,24 +470,33 @@ void Update::reset_timestep(int narg, char **arg)
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
reset timestep reset timestep
called from rerun command and input script (indirectly) called from input script (indirectly) or rerun command
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Update::reset_timestep(bigint newstep) void Update::reset_timestep(bigint newstep)
{ {
if (newstep < 0) error->all(FLERR,"Timestep must be >= 0");
bigint oldstep = ntimestep;
ntimestep = newstep; ntimestep = newstep;
if (ntimestep < 0) error->all(FLERR,"Timestep must be >= 0");
// set atimestep to new timestep // if newstep >= oldstep, update simulation time accordingly
// so future update_time() calls will be correct // if newstep < oldstep, zero simulation time
atimestep = ntimestep; if (newstep >= oldstep) update_time();
// trigger reset of timestep for output if (newstep < oldstep) {
// do not allow any timestep-dependent fixes to be already defined atime = 0.0;
atimestep = newstep;
}
// changes to output that depend on timestep
// no active dumps allowed
output->reset_timestep(ntimestep); output->reset_timestep(ntimestep);
// do not allow timestep-dependent fixes to be defined
for (const auto &ifix : modify->get_fix_list()) for (const auto &ifix : modify->get_fix_list())
if (ifix->time_depend) if (ifix->time_depend)
error->all(FLERR, "Cannot reset timestep with time-dependent fix {} defined",ifix->style); error->all(FLERR, "Cannot reset timestep with time-dependent fix {} defined",ifix->style);
@ -508,7 +517,7 @@ void Update::reset_timestep(bigint newstep)
if (icompute->timeflag) icompute->clearstep(); if (icompute->timeflag) icompute->clearstep();
} }
// Neighbor Bin/Stencil/Pair classes store timestamps that need to be cleared // neighbor Bin/Stencil/Pair classes store timestamps that need to be cleared
neighbor->reset_timestep(ntimestep); neighbor->reset_timestep(ntimestep);
} }

9
unittest/formats/test_dump_atom.cpp
View File

@ -74,6 +74,13 @@ public:
END_HIDE_OUTPUT(); END_HIDE_OUTPUT();
} }
void close_dump()
{
BEGIN_HIDE_OUTPUT();
command("undump id");
END_HIDE_OUTPUT();
}
void generate_text_and_binary_dump(std::string text_file, std::string binary_file, void generate_text_and_binary_dump(std::string text_file, std::string binary_file,
std::string dump_modify_options, int ntimesteps) std::string dump_modify_options, int ntimesteps)
{ {
@ -505,6 +512,7 @@ TEST_F(DumpAtomTest, rerun)
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
ASSERT_EQ(count_lines(dump_file), 82); ASSERT_EQ(count_lines(dump_file), 82);
continue_dump(1); continue_dump(1);
close_dump();
lmp->output->thermo->evaluate_keyword("pe", &pe_2); lmp->output->thermo->evaluate_keyword("pe", &pe_2);
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
ASSERT_EQ(count_lines(dump_file), 123); ASSERT_EQ(count_lines(dump_file), 123);
@ -532,6 +540,7 @@ TEST_F(DumpAtomTest, rerun_bin)
lmp->output->thermo->evaluate_keyword("pe", &pe_1); lmp->output->thermo->evaluate_keyword("pe", &pe_1);
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
continue_dump(1); continue_dump(1);
close_dump();
lmp->output->thermo->evaluate_keyword("pe", &pe_2); lmp->output->thermo->evaluate_keyword("pe", &pe_2);
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
HIDE_OUTPUT([&] { HIDE_OUTPUT([&] {

10
unittest/formats/test_dump_custom.cpp
View File

@ -73,6 +73,13 @@ public:
END_HIDE_OUTPUT(); END_HIDE_OUTPUT();
} }
void close_dump()
{
BEGIN_HIDE_OUTPUT();
command("undump id");
END_HIDE_OUTPUT();
}
void generate_text_and_binary_dump(std::string text_file, std::string binary_file, void generate_text_and_binary_dump(std::string text_file, std::string binary_file,
std::string fields, std::string dump_modify_options, std::string fields, std::string dump_modify_options,
int ntimesteps) int ntimesteps)
@ -330,6 +337,7 @@ TEST_F(DumpCustomTest, rerun)
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
ASSERT_EQ(count_lines(dump_file), 82); ASSERT_EQ(count_lines(dump_file), 82);
continue_dump(1); continue_dump(1);
close_dump();
lmp->output->thermo->evaluate_keyword("pe", &pe_2); lmp->output->thermo->evaluate_keyword("pe", &pe_2);
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
ASSERT_EQ(count_lines(dump_file), 123); ASSERT_EQ(count_lines(dump_file), 123);
@ -338,6 +346,7 @@ TEST_F(DumpCustomTest, rerun)
}); });
lmp->output->thermo->evaluate_keyword("pe", &pe_rerun); lmp->output->thermo->evaluate_keyword("pe", &pe_rerun);
ASSERT_DOUBLE_EQ(pe_1, pe_rerun); ASSERT_DOUBLE_EQ(pe_1, pe_rerun);
HIDE_OUTPUT([&] { HIDE_OUTPUT([&] {
command(fmt::format("rerun {} first 2 last 2 every 1 post yes dump x y z", dump_file)); command(fmt::format("rerun {} first 2 last 2 every 1 post yes dump x y z", dump_file));
}); });
@ -359,6 +368,7 @@ TEST_F(DumpCustomTest, rerun_bin)
lmp->output->thermo->evaluate_keyword("pe", &pe_1); lmp->output->thermo->evaluate_keyword("pe", &pe_1);
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
continue_dump(1); continue_dump(1);
close_dump();
lmp->output->thermo->evaluate_keyword("pe", &pe_2); lmp->output->thermo->evaluate_keyword("pe", &pe_2);
ASSERT_FILE_EXISTS(dump_file); ASSERT_FILE_EXISTS(dump_file);
HIDE_OUTPUT([&] { HIDE_OUTPUT([&] {